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ABSTRACT 


Combining  horizontal  flow  treatment  wells  (HFTWs)  with  in  situ  biodegradation 
is  an  innovative  approach  with  the  potential  to  remediate  perchlorate-contaminated 
groundwater.  A  technology  model  was  recently  developed  that  combines  the 
groundwater  flow  induced  by  HFTWs  with  in  situ  biodegradation  processes  that  result 
from  using  the  HFTWs  to  mix  electron  donor  into  perchlorate-contaminated 
groundwater.  A  field  demonstration  of  this  approach  is  planned  to  begin  this  year. 

In  order  to  apply  the  technology  in  the  field,  project  managers  need  to  understand 
how  contaminated  site  conditions  and  technology  design  parameters  impact  technology 
performance.  One  way  to  gain  this  understanding  is  to  use  the  technology  model  to  select 
engineering  design  parameters  that  optimize  performance  under  given  site  conditions.  In 
particular,  a  project  manager  desires  to  design  a  system  that  1)  maximizes  perchlorate 
destruction,  2)  minimizes  treatment  expense,  and  3)  attains  regulatory  limits  on 
downgradient  contaminant  concentrations.  Unfortunately,  for  a  relatively  complex 
technology  like  in  situ  bioremediation,  with  a  number  of  engineering  design  parameters 
to  determine,  as  well  as  multiple  objectives,  system  optimization  is  not  straightforward. 

In  this  study,  a  multi-objective  genetic  algorithm  (MOGA)  is  used  to  determine 
design  parameter  values  (flow  rate,  well  spacing,  concentration  of  injected  electron 
donor,  and  injection  schedule)  that  optimize  the  first  two  objectives  noted;  to  maximize 
perchlorate  destruction  while  minimizing  cost.  Four  optimization  runs  are  performed, 
using  two  different  remediation  time  spans  (300  and  600  days)  for  two  different  sets  of 
site  conditions.  Results  from  all  four  optimization  runs  indicate  that  the  relationship 


IX 


between  perchlorate  mass  removal  and  operating  cost  is  positively  correlated  and 
nonlinear.  For  equivalent  operating  times  and  costs,  the  optimized  solutions  show  that,  as 
expected,  the  technology  achieves  higher  mass  removals  for  the  site  having  both  higher 
hydraulic  conductivity  and  higher  initial  source  concentration.  Results  from  all  four  runs 
show  that  increased  perchlorate  mass  removal  is  not  necessarily  correlated  with 
diminished  downgradient  perchlorate  concentration,  suggesting  that  it  may  be  important 
to  incorporate  minimization  of  downgradient  perchlorate  concentration  as  an  additional 
objective  or  constraint  in  the  multi-objective  optimization  scheme. 

The  optimization  software  developed  in  this  study  can  serve  as  a  tool  for  both 
optimizing  future  applications  of  this  innovative  bioremediation  technology,  and  helping 
us  to  better  understand  how  HFTWs  can  be  used  in  conjunction  with  in  situ 
biodegradation.  This  study  contributes  to  efforts  taken  to  resolve  groundwater 
contamination  problems  caused  by  perchlorate  releases  across  the  United  States. 
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OPTIMIZING  AN  IN  SITU  BIOREMEDIATION  TECHNOLOGY  TO  MANAGE 


PERCHLORATE-CONTAMINATED  GROUNDWATER 


1.0  INTRODUCTION 


1.1  BACKGROUND 

Perchlorate  is  an  oxyanion  that  the  aerospace  industry  has  used  since  the  1940’s 
as  a  constituent  in  solid  rocket  fuel  (EPA,  1999:1;  Herman  and  Frankenberger, 
1998:750).  Due  to  the  absence  of  legal  restrictions,  a  lack  of  knowledge  of  perchlorate 
health  effects,  and  a  deficient  understanding  of  the  processes  affecting  perchlorate  fate 
and  transport,  high  levels  of  ammonium  perchlorate  were  discharged  into  the 
environment  (Urbansky,  1998:82),  resulting  in  perchlorate  contamination  problems  that 
we  face  today.  Perchlorate  contamination  from  these  past  practices  now  affects  the 
drinking  water  of  15  million  U.S.  citizens  (EPA,  1999:1);  this  contamination  is 
particularly  significant  in  western  states  like  California,  Utah,  and  Arizona  (Urbansky, 
1998:82). 

The  chief  health  problem  caused  by  perchlorate  is  due  to  its  potential  to  interfere 
with  hormone  production  in  humans.  The  thyroid  gland  normally  uptakes  iodide  from 
the  bloodstream  to  make  hormones;  however,  the  presence  of  perchlorate  in  the 
bloodstream  causes  the  thyroid  gland  to  uptake  perchlorate  instead  of  iodide,  thereby 
disrupting  hormone  production.  Animal  studies  also  show  perchlorate's  potential  to 
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interfere  with  muscle  movement  (Urbansky,  1998:83).  The  presence  of  perchlorate  in  the 
environment  triggers  even  more  concern  because  the  consequences  of  prolonged,  low- 
dose  exposure  have  yet  to  be  ascertained  (Logan,  2001:484A). 

Unfortunately,  perchlorate  is  mobile  and  persistent  in  the  natural  environment. 
According  to  Flowers  and  Hunt  (2000:177),  perchlorate  is  "expected  to  be  highly  mobile 
in  surface  and  groundwaters".  Under  nearly  neutral  pH  conditions,  which  are  found  in 
most  groundwaters,  perchlorate  would  probably  not  sorb  to  mineral  surfaces  (Flowers 
and  Hunt,  2000:177);  additionally,  perchlorate  ions  are  "not  retarded  during  groundwater 
transport"  (Logan,  2001 :483A).  Natural  destruction  of  perchlorate  is  unlikely.  Dissolved 
perchlorate  salts  not  only  resist  reaction  via  coordination,  but  they  also  "hardly  react  at  all 
in  any  manner" (Espenson,  2000: 1).  Furthermore,  insufficient  proof  exists  that 
perchlorate  would  naturally  degrade  via  biological  transfonnation  (Flowers  and  Hunt, 
2000:177). 

The  perchlorate  problem  is  exacerbated  because  remediation  of  perchlorate- 
contaminated  water  is  difficult.  Logan  (2001 :484A)  asserts  that  there  is  “no  obvious 
treatment  technology  for  removing  perchlorate  from  water”.  Ion  exchange,  air  stripping, 
carbon  adsorption,  and  advanced  oxidation  do  not  provide  cost-efficient  perfonnance 
(Logan,  1998:70).  Physical  removal  methods  do  not  destroy  perchlorate;  they  merely 
concentrate  perchlorate  elsewhere,  which  creates  waste  disposal  problems  (Urbansky, 
1998:90).  Pump-and-treat  systems  are  poor  remediation  choices  because  they  necessitate 
pumping  contaminated  water  to  the  surface,  which  increases  treatment  costs  and 
introduces  risk  of  exposure  to  perchlorate  (Ferland  and  Goltz,  2001:45). 
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In  situ  biodegradation,  however,  offers  a  potential  solution  to  the  perchlorate 
problem.  According  to  Logan  (2001:486A  -  487A),  injection  of  substrates  directly  into 
the  subsurface  can  promote  in  situ  microbial  degradation;  such  a  strategy  was  successful 
in  perchlorate  remediation  projects  at  MacGregor,  TX,  and  Sacramento,  CA. 

A  critical  issue  in  injecting  chemicals  to  promote  in  situ  biodegradation  is 
ensuring  the  injected  chemical  sufficiently  mixes  with  the  contaminated  groundwater.  A 
new  technology  that  uses  so-called  horizontal  flow  treatment  wells  (HFTWs)  to  effect 
mixing  to  promote  in  situ  biodegradation  has  potential  to  remediate  perchlorate 
contamination  (McCarty  et  al,  1998;  Parr,  2002).  A  field  test  at  Edwards  AFB,  CA, 
successfully  employed  a  similar  technology  to  remediate  TCE-contaminated 
groundwater;  in  that  test,  pulses  of  toluene,  oxygen,  and  hydrogen  peroxide  were  injected 
into  groundwater  that  was  circulating  between  two  HFTWs  (McCarty  et  al.,  1998:88). 
The  Environmental  Security  Technology  Certification  Program  (ESTCP)  funded  a 
project  to  demonstrate  that  injection  and  in  situ  mixing  of  electron  donor  into  perchlorate- 
contaminated  groundwater  using  HFTWs  could  be  an  effective  treatment  technology 
(ESTCP,  2002).  Parr  (2002)  developed  a  model  that  combines  HFTW-induced 
groundwater  flow  with  in  situ  biodegradation  processes  that  result  from  use  of  HFTWs  to 
mix  electron  donor  into  perchlorate-contaminated  groundwater. 

To  apply  this  new  technology  in  the  field,  project  managers  must  understand  how 
contaminated  site  conditions  and  design  parameters  affect  technology  performance.  An 
approach  to  gaining  this  understanding  is  to  use  a  technology  model  in  order  to  select 
engineering  design  parameters  that  optimize  performance  for  given  site  conditions.  In 
particular,  a  project  manager  desires  a  system  design  that  1)  maximizes  perchlorate 
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destruction,  2)  minimizes  treatment  expense,  and  3)  attains  downgradient  regulatory 
requirements. 

The  dilemma  of  pursuing  separate  objectives  results  in  a  multi-objective 
optimization  problem.  In  multi-objective  optimization,  the  decision  maker  must 
manipulate  decision  variables  ( i.e .  engineering  parameters)  to  ideally  obtain  a  solution 
that  optimizes  all  objectives.  Furthennore,  this  optimum  solution  must  be  feasible  (i.e. 
satisfy  constraints  relevant  to  the  problem). 

Unfortunately,  a  single  feasible  solution  that  optimizes  all  objectives  usually 
doesn’t  exist.  In  fact,  most  multi-objective  problems  have  inherent  trade-off  dilemmas, 
where  improving  perfonnance  in  one  objective  usually  entails  worsening  performance  in 
another.  For  example,  the  decision  maker  may  discover  that  reducing  the  operating  cost 
of  the  remediation  technology  (desirable)  incurs  degraded  performance  in  contaminant 
destruction  (undesirable).  When  tradeoffs  exist  among  competing  objectives,  the  solution 
to  a  multi-objective  problem  involves  finding  multiple  solutions  that  are  nondominated, 
or  Pareto  optimal.  If  a  solution  x  is  Pareto  optimal,  then  no  other  solutions  exist  that 
perform  better  than  x  across  all  objectives.  Pareto  optimal  solutions  comprise  the  Pareto 
optimal  set  (P*),  and  they  exhibit  a  tradeoff  relationship  in  objective  space  called  the 
Pareto  front  (PF*)  (Coello  Coello  et  al.,  2002:10-12). 

Simultaneously  optimizing  perfonnance  and  cost  for  this  innovative  in  situ 
perchlorate  bioremediation  technology  is  a  complicated  undertaking;  intense 
computational  effort  is  necessary  to  handle  the  nonlinearities  of  the  model’s  equations 
and  account  for  constraints  while  choosing  the  best  combination  of  several  design 
variables,  each  of  which  can  vary  over  a  wide  range.  Hence,  in  order  to  optimize 
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application  of  the  technology  and  gain  greater  understanding  of  technology  costs  and 
performance,  it  becomes  necessary  to  select  an  optimization  technique  that  can  handle 
multiple  objectives.  One  viable  technique  is  the  genetic  algorithm  (GA),  a  stochastic 
search  algorithm  whose  operators  mimic  biological  genetics.  GAs  have  the  ability  to 
operate  on  multiple  candidate  solutions,  as  opposed  to  only  one  solution,  which  renders 
them  conducive  to  multi-objective  optimization. 


1.2  RESEARCH  GOALS  AND  OBJECTIVES 

The  overall  goal  of  this  thesis  effort  is  to  enhance  understanding  of  how  in  situ 
bioremediation,  used  in  conjunction  with  HFTWs,  can  e  applied  to  manage  perchlorate- 
contaminated  groundwater.  This  goal  can  be  divided  into  the  following  specific  research 
objectives: 

•  Develop  a  multi-objective  genetic  algorithm  (MOGA)  that  can  be  used  in 
conjunction  with  the  in  situ  perchlorate  bioremediation  technology  model 

•  Determine  values  of  design  variables  that  optimize  technology  cost  and 
performance  given  different  time  periods  and  contaminated-site  conditions 

•  Assess  how  the  technology  perfonns  in  reducing  downgradient  concentration 
given  different  time  periods  and  contaminated-site  conditions 
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1.3  RESEARCH  APPROACH 


•  Formulate  a  multi-objective  problem  (objective  functions,  decision  variables, 
constraints)  for  application  of  Parr’s  (2002)  technology  model  to  manage 
perchlorate-contaminated  groundwater. 

•  Design  and  encode  an  appropriate  MOGA  based  on  a  review  of  the  literature 

•  Couple  the  technology  model  with  the  developed  MOGA 

•  Determine  optimal  technology  design  parameters  for  given  site  conditions 


1.4  SCOPE  AND  LIMITATIONS  OF  RESEARCH 

•  This  study  assumes  the  technology  model  developed  by  Parr  (2002)  adequately 
simulates  remediation  of  perchlorate-contaminated  groundwater  through  electron 
donor  injection  by  HFTWs. 

•  Optimization  scenarios  only  involve  a  two-well  HFTW  system  to  facilitate 
computation.  It  is  assumed  results  obtained  from  the  two-well  system  can  be 
scaled  up  to  manage  wider  contaminant  plumes. 

•  During  optimization,  engineering  decision  variables  are  restricted  to  a  range  of 
values  appropriate  for  the  formulated  problem. 

•  MOGAs  use  probabilistic,  rather  than  deterministic,  operators  to  explore  the 
search  space  for  Pareto-optimal  solutions.  Hence,  both  the  quantity  and  diversity 
of  Pareto  optimal  points  that  a  MOGA  can  find  are  inherently  uncertain. 

•  This  study  limits  the  multi-objective  problem  to  two  objectives  for  the  technology 
model:  1)  mass  removal  of  contaminant  and  2)  operating  cost. 
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This  study  assumes  that  the  operating  cost  for  the  technology  depends  solely  upon 
pump  operation  and  injection  of  electron  donor. 


2.0  LITERATURE  REVIEW 


2.1  INTRODUCTION 

This  chapter  begins  by  presenting  information  on  perchlorate,  including  health  effects, 
regulatory  issues,  contaminant  characteristics,  and  remediation  strategies.  The  next 
section  discusses  horizontal  flow  treatment  wells  (HFTWs),  which  may  offer  a  means  of 
remediating  perchlorate-contaminated  sites.  Section  2.4  introduces  the  reader  to  multi¬ 
objective  optimization,  with  emphasis  on  Pareto  optimality.  Section  2.5  provides  an 
overview  of  genetic  algorithms  (GAs).  Finally,  Section  2.6  addresses  how  GAs  may  be 
applied  to  solve  multi-objective  problems. 


2.2  PERCHLORATE 

Perchlorate  (CKV)  salts  are  used  by  various  industries,  including  the  auto, 
chemical,  rubber,  and  fabric  industries.  However,  the  most  notable  user  is  the  aerospace 
industry.  In  the  mid- 1 940 ’s,  the  aerospace  industry  started  large-scale  production  of 
ammonium  perchlorate  (NH4CIO4)  as  an  oxidizer  for  solid  rocket  fuel  (EPA,  1999:1). 
The  relatively  short  shelf  life  of  this  chemical  and  absence  of  disposal  regulations 
resulted  in  recurring  discharges  of  ammonium  perchlorate  to  the  environment  (EPA, 
1999:1;  Urbansky,  1998:82).  More  recently,  demilitarization  has  led  to  missile 
disassembly,  which  may  also  account  for  discharges  (Urbansky,  1998:82). 

Solid  salts  containing  perchlorate  dissolve  readily  in  water.  In  particular, 
ammonium  perchlorate  “is  highly  soluble  and  dissociates  completely”  to  ammonium 
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(NH4+)  and  perchlorate  (CIO4')  ions  (Urbansky,  1998:82).  Free  perchlorate  ions  result  in 
perchlorate  contamination,  which  is  now  a  pervasive  problem  in  the  United  States. 
Eighteen  states  have  confirmed  contamination  (Logan,  200 1:483 A);  California,  Nevada, 
and  Arizona  alone  have  contaminated  water  supplies  impacting  15  million  citizens  (EPA, 
1999:1). 

Perchlorate  is  problematic  because  of  its  potential  to  interfere  with  hormone 
production  in  humans.  Under  nonnal  circumstances,  the  thyroid  gland  extracts  iodide 
from  the  bloodstream  to  produce  hormones.  Because  the  thyroid  gland  has  higher 
selectivity  for  perchlorate  than  iodide,  however,  perchlorate  is  likely  to  displace  iodide  in 
the  gland,  thereby  disrupting  hormone  production  (Urbansky,  1998:83).  The  U.S.  EPA 
(Environmental  Protection  Agency)  observed  health  effects  of  perchlorate  that  was 
chemotherapeutically  administered  to  patients  afflicted  with  hyperthyroidism;  after  2 
months,  minimum  doses  of  6  mg/kg/day  caused  “fatal  bone  marrow  changes”  (Urbansky, 
1998:83).  Also,  animal  studies  show  that  high  perchlorate  concentrations  can  disrupt 
muscle  movement  (Urbansky,  1998:83).  Uncertainty  exists  as  to  toxicity  of  perchlorate, 
particularly  due  to  prolonged,  low-dose  exposure  (EPA,  1999:1;  Logan,  2001:484A). 

Until  1997,  perchlorate  concentrations  below  100  parts  per  billion  (ppb)  were 
undetectable.  In  April  1997,  the  California  Department  of  Health  Services  (CDHS) 
developed  a  method  to  detect  perchlorate  levels  as  low  as  4  ppb.  Subsequent  to  the 
development  of  this  new  capability,  perchlorate  has  been  detected  in  the  following  states: 
AR,  AZ,  CA,  IA,  IN,  KS,  MD,  NM,  NV,  NY,  PA,  TX,  UT,  and  WV  (EPA,  1999: 1).  In 
1997  both  California  and  Nevada  set  action  levels  for  perchlorate  at  18  ppb;  in  1999 
Arizona  and  Texas  set  action  levels  of  31  and  22  ppb,  respectively  (EPA,  1999:2). 
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Although  the  United  States  Environmental  Protection  Agency  (USEPA)  placed 
perchlorate  on  the  1998  Contaminant  Candidate  List,  USEPA  still  has  not  issued  federal 
drinking  water  regulations  for  perchlorate.  However,  on  January  2,  2002,  EPA  did 
release  a  revised  draft  of  its  toxicity  assessment  for  public  review.  The  draft  assessment, 
entitled  "Perchlorate  Environmental  Contamination:  Toxicological  Review  and  Risk 
Characterization,"  is  a  culmination  of  research  efforts  since  1997  pertaining  to 
perchlorate’s  health  impacts  (Federal  Register,  2002:75-76).  As  of  this  date,  however, 
the  draft  assessment  is  still  in  review. 

The  California  Environmental  Protection  Agency  (CEP A)  proposed  a  public 
health  goal  (PHG)  of  6  ppb  perchlorate  in  drinking  water  based  on  a  No  Observed 
Adverse  Effect  Level  (NOAEL)  of  10  pg/kg-day  and  uncertainty  factor  of  30  (CEP A, 
2002:2). 

Perchlorate  anions  persist  in  the  natural  environment.  A  perchlorate  anion  (CIO4') 
consists  of  one  chlorine  atom  centered  among  four  oxygen  atoms  in  a  tetrahedral 
structure,  with  the  negative  charge  evenly  distributed  among  the  oxygen  atoms 
(Espenson,  2000:1).  This  symmetric  charge  distribution  makes  the  perchlorate  ion 
resistant  to  complexation  with  metals  (Espenson,  2000: 1).  Thus,  even  though  perchlorate 
is  a  strong  oxidant,  “(p)erchlorate  reactions  demonstrate  high  kinetic  barriers”  that  render 
the  ion  virtually  unreactive  and,  therefore,  persistent  (Espenson,  2000:2).  There  appears 
to  be  insufficient  evidence  that  perchlorate  is  biologically  transformable  “under  natural 
conditions”  (Flowers  and  Hunt,  2000:177). 

In  addition  to  persistence,  perchlorate  also  demonstrates  excellent  mobility  in 
groundwater.  Reported  solubility  of  ammonium  perchlorate  is  relatively  high  at  200  g/L, 
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and  “sodium,  calcium,  and  magnesium  salts  are  even  more  soluble”  (Flowers  and  Hunt, 
2000:177).  Under  nearly  neutral  pH  conditions  typical  of  most  groundwater,  sorption  of 
perchlorate  to  mineral  surfaces  is  unlikely  (Flowers  and  Hunt,  2000:177)  and  perchlorate 
ions  are  "not  retarded  during  groundwater  transport"  (Logan,  200 1:483 A). 

Methods  to  treat  perchlorate-contaminated  groundwater  can  be  implemented  in 
two  ways:  ex  situ  (above  ground)  or  in  situ  (in  place,  below  ground).  Ex  situ  treatment 
necessitates  application  of  a  pump-and-treat  (PAT)  system  to  pump  the  contaminated 
water  from  the  subsurface  to  the  surface  where  treatment  occurs.  Disadvantages  of  PAT 
are  discussed  in  Section  2.3.  However,  ex  situ  treatment  is  advantageous  because  it 
allows  implementation  of  various  conventional  and  innovative  treatment  methods,  such 
as  anion  exchange,  reverse  osmosis  (RO),  chemical/electrochemical  reduction,  and 
engineered  biotreatment.  A  more  detailed  discussion  of  each  of  these  methods  follows. 

Anion  exchange  is  an  existing  technology  that  uses  a  resin  to  extract  aqueous 
perchlorate  and  replace  it  with  an  innocuous  anion  (Urbansky  and  Schock,  1999:86). 
Unfortunately,  current  resins  have  poor  selectivity  for  perchlorate,  and  resins  with 
enhanced  selectivity  are  expensive  (Urbansky  and  Schock,  1999:84).  Selectivity  for 
perchlorate  worsens  at  lower  perchlorate  concentrations  because  resins  tend  to  favor 
more  abundant  anions  for  exchange  (Urbansky  and  Schock,  1999:84).  Resins  can 
remove  aqueous  perchlorate  only  until  an  equilibrium  concentration  occurs,  after  which 
further  removal  is  impossible  (Urbansky  and  Schock,  1999:84).  Time  is  lost  regenerating 
the  resin  (Urbansky  and  Schock,  1999:84),  and  brines  used  to  regenerate  resins  can  have 
perchlorate  concentrations  high  enough  to  warrant  disposal  concerns  (Logan, 

200L484A). 
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In  RO,  high  pressure  forces  untreated  water  through  a  semiporous  polymer 
membrane  (Urbansky  and  Schock,  1999:86).  The  membrane,  which  is  impermeable  to 
dissolved  salts,  acts  as  a  fdter  (Urbansky  and  Schock,  1999:86).  Water  passes  through 
the  membrane  and  dissolved  salts  stay  behind,  so  that  the  membrane-filtered  water  is 
relatively  deionized  (Urbansky  and  Schock,  1999:87).  However,  RO  membranes  are 
susceptible  to  fouling  from  metallic  compounds,  natural  organic  matter,  and  microbes; 
such  fouling  necessitates  costly  membrane  replacement  (Urbansky,  1998:89;  Urbansky 
and  Schock,  1999:86-87). 

Because  anion  exchange  and  RO  are  both  physical  removal  methods,  they  have  a 
common  drawback:  waste  disposal  (Urbansky  and  Schock,  1999:85;  Urbansky,  1998:90). 
Physical  removal  merely  transfers  perchlorate  from  water  to  another  medium  instead  of 
destroying  it.  The  transferred  perchlorate  must  either  undergo  further  treatment  or 
disposal.  As  Urbansky  (1998)  says,  “Although  these  techniques  take  the  perchlorate  out, 
they  concentrate  it  somewhere  else  where  it  must  be  dealt  with  later”  (90). 

Ex  situ  perchlorate  destruction  is  possible  with  innovative  technologies,  such  as 
electrochemical  reduction,  reduction  using  titanous  ions  or  metallic  iron/UV  light,  and 
biological  treatment.  In  electrochemical  reduction,  a  cathode  applies  electric  current 
directly  to  water  at  a  high  potential  to  reduce  perchlorate  to  chloride  (Urbansky  and 
Schock,  1999:85).  However,  movement  of  perchlorate  ions  toward  the  electrode  can  take 
a  long  time,  and  electrodes  are  susceptible  to  corrosion  and  fouling  from  natural  organic 
matter  (Urbansky  and  Schock,  1999:85). 

Earley  et  al.  (2000)  hypothesized  that  perchlorate  might  be  effectively  destroyed 
by  reaction  with  trivalent  titanous  ions  [Ti(H20)6  ]  in  an  ethanol  media,  which  catalyze 
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the  reaction.  Amadei  and  Earley  (2001)  also  proposed  catalysts  that  offer  faster  reaction 
rates  than  ethanol  media.  Parr  (2002)  points  out,  however,  that  this  technology  has  “very 
limited  laboratory  data,”  and  pilot-scale  tests  have  yet  to  be  accomplished  (63). 

According  to  Gurol  and  Kim  (2000),  exposing  perchlorate  to  UV  light  and 
metallic  iron  (Fe°)  under  anoxic  conditions  can  transform  perchlorate  to  chloride  and 
water.  It  is  hypothesized  that  metallic  iron  first  adsorbs  the  perchlorate  ion  and  then 
undergoes  oxidation,  with  UV  light  expediting  the  electron  transfer  from  the  iron  to  the 
perchlorate.  However,  high  reaction  rates  were  achievable  only  with  very  high  UV 
intensity  (Gurol  and  Kim,  2000).  Limited  data  on  the  metallic  iron/UV  technology 
suggest  this  technology  is  currently  inappropriate  for  lowering  perchlorate  levels  to 
below  regulatory  limits  (Parr,  2002:63). 

Research  on  ex  situ  engineered  biological  treatment  (“biotreatmenf ’)  of 
perchlorate-contaminated  water  has  largely  focused  on  flowing  water  through  a  column- 
type  bioreactor  to  effect  treatment.  In  general,  these  bioreactors  contain  microbial 
support  media  (GAC,  sand,  plastic)  that  are  inoculated  with  perchlorate-reducing 
microbes  (PRMs).  Inflow  is  a  mixture  of  perchlorate  and  substrate  (electron  donor). 
Inside  the  reactor,  PRMs  reduce  perchlorate  via  the  following  pathway  (Rikken  el  al., 
1996:425): 

ClOfi  (perchlorate)  — >  CIO3"  (chlorate)  — >  CIO2"  (chlorite)  — >  Cl"  (chloride)  +  O2 
In  this  sequential  reduction,  the  substrate  acts  as  an  electron  donor,  and  perchlorate  acts 
as  the  terminal  electron  acceptor. 

The  two  primary  designs  of  ex  situ  bioreactors  are  fluidized-  and  fixed-bed. 
Fluidized-bed  reactors  use  high  flow  rates  to  mix  the  support  medium  (GAC,  sand) 
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within  the  reactor  (Logan,  200 1:485 A).  In  fixed-bed  reactors,  the  support  medium 
remains  stationary  (“fixed”)  during  water  treatment  (Logan,  200 1:485 A).  Pilot  versions 
of  fluidized-  and  fixed-bed  bioreactors  have  demonstrated  the  ability  to  reduce 
perchlorate  concentrations  to  below  the  current  detection  limit  (4  ppb).  Full-scale 
fluidized-bed  reactors  in  California  have  effectively  reduced  perchlorate  levels  to  <  4  ppb 
as  well  (Hatzinger  et  al.,  2000:7;  Greene  and  Pitre,  2000:252). 

In  situ  methods  offer  the  ability  to  treat  perchlorate-contaminated  water  while  it  is 
still  in  the  ground.  In  situ  methods  are  favorable  because,  unlike  ex  situ  methods,  they  do 
not  require  pumping  groundwater  to  the  surface  for  treatment,  thereby  reducing  both 
operating  expenses  and  human  exposure  to  the  contaminant.  However,  in  situ  treatment 
technologies  are  limited  to  those  technologies  that  can  be  applied  in-well  or  below 
ground.  Parr  (2002)  surveyed  perchlorate  treatment  technologies  and  determined  that 
biotreatment  is  the  only  method  readily  applicable  for  in  situ  use. 

Both  in  situ  and  ex  situ  biotreatment  rely  on  the  same  reaction  pathway  (Rikken  et 
al.,  1996:425)  to  desynthesize  perchlorate.  However,  in  situ  biotreatment  requires  1) 
presence  of  PRMs  at  the  contaminated  site  and  2)  a  means  of  creating  anoxic  conditions 
for  perchlorate  reduction  (Logan,  2001 :486A).  The  first  requirement  is  perhaps  the 
easiest  to  fulfill.  PRMs  appear  to  be  found  in  soil  at  perchlorate-contaminated  sites  and 
are  “widely  distributed  in  nature”  (Wu  et  al.,  2001 : 1 19,  125).  To  fulfill  the  second 
requirement,  Logan  (2001)  proposes  either  use  of  biobarriers  or  direct  injection  of 
substrates  into  the  ground  (486A-487A). 

Biobarriers  may  be  a  cost-effective  way  of  implementing  in  situ  biotreatment 
(McMaster  et  al.,  2001 :301).  A  biobarrier  contains  “high  concentrations  of  organic 
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matter”  and  establishes  a  vertical  zone  of  bioactivity  in  the  ground  (Logan,  200L486A- 
487A;  Domenico  and  Schwarz,  1998:450).  As  contaminated  groundwater  passes  through 
the  barrier,  microorganisms  chemically  transfonn  the  target  contaminant  (Domenico  and 
Schwarz,  1998:450).  A  biobarrier  is  a  “passive”  remediation  technology  that  can  remain 
in-ground  for  years  with  negligible  maintenance  (Domenico  and  Schwarz,  1998:450).  In 
MacGregor,  TX,  biobarriers  were  used  to  remediate  soil  contaminated  with  high  levels  of 
perchlorate  (Logan,  200L487A).  Biobarriers  were  installed  by  digging  trenches  and 
filling  them  with  organic  materials  and  gravel  (Logan,  2001 :487A).  By  directing  water 
flow  through  the  biobarrier,  perchlorate  dropped  from  27,000  pg/L  to  nondetectable 
levels  (Logan,  200L487A).  As  a  passive  technology,  however,  biobarriers  may  be 
bypassed  as  groundwater  flow  conditions  change  over  time.  In  addition,  biobarriers  can 
only  be  emplaced  to  a  limited  depth. 

Injection  of  substrates  into  the  ground  is  another  in  situ  biotreatment  strategy 
(Logan,  2001 :486A-487A).  In  May  2000  a  field  demonstration  at  the  Aerojet  Superfund 
Site  in  Sacramento,  CA,  showed  the  effectiveness  of  in  situ  biodegradation  of  perchlorate 
(McMaster  et  al.,  200 1 :297).  The  test  site  had  a  perchlorate  plume  with  concentrations 
ranging  from  10,000  -  15,000  pg/L.  A  closed-loop  recirculation  system,  located  within 
the  plume’s  center,  extracted  water  from  the  aquifer,  added  time-pulsed  doses  of  acetate, 
and  injected  the  mixture  back  into  the  aquifer  to  stimulate  bioactivity.  The  system 
reduced  perchlorate  levels  to  <18  pg/L,  the  CDHS  action  level  for  perchlorate  in  drinking 
water  (McMaster  et  al.,  2001 :297,  299).  Note  however,  that  the  closed-loop  recirculation 
system  applied  at  the  Aerojet  Site  required  extraction  and  reinjection  of  groundwater, 
using  a  strategy  similar  to  PAT.  In  order  to  effectively  achieve  injection  and  mixing  of 
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substrate  into  perchlorate-contaminated  groundwater,  and  delivery  of  the  mixture  to 
indigenous  PRMs,  without  the  need  to  pump  water  to  the  surface,  an  innovative 
technology  has  been  proposed  (Parr,  2002).  This  technology,  horizontal  flow  treatment 
wells  (HFTWs),  is  discussed  in  the  next  section. 


2.3  HORIZONTAL  FLOW  TREATMENT  WELL  (HFTW)  SYSTEM 

The  operating  concept  for  an  HFTW  system  is  for  the  treatment  wells  to  be 
installed  at  the  distal  end  of  a  contaminant  plume.  The  plume  would  be  captured  by  the 
wells  and  the  contaminated  groundwater  treated  (i.e.  contaminant  mass  destroyed)  so  that 
contaminant  concentrations  downgradient  of  the  treatment  wells  would  meet  regulatory 
requirements  (Ferland  and  Goltz,  2001:46).  Figures  1  and  2  show  cross-sectional  and 
plan  views,  respectively,  of  HFTW  operation. 

As  shown  in  Figure  1,  an  HFTW  has  two  screens,  each  located  in  a  different 
subsurface  horizon.  The  HFTW  extracts  contaminated  groundwater  through  a  screen  in 
one  horizon,  treats  the  water,  and  then  discharges  the  water  through  the  screen  in  the 
other  horizon.  In-well  treatment  can  involve  the  application  of  a  physical,  chemical,  or 
biological  process  to  remove  or  destroy  the  contaminant  (Gandhi  et  al.,  2002a:4);  or,  as 
depicted  in  the  figure,  it  can  involve  injection  of  chemicals  into  the  flowing  contaminated 
groundwater  in  order  to  establish  bioactive  zones  outside  the  wells’  discharge  screens, 
where  the  contaminant  is  biodegraded  by  indigenous  microorganisms  (McCarty  et  al., 
1998:90).  HFTWs  are  installed  in  pairs,  with  one  well  pumping  in  an  “up flow  mode” 
while  the  adjacent  well  pumps  in  a  “downflow  mode,”  resulting  in  water  recirculating 
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between  the  two  wells  in  a  horizontal  flow  pattern  (Ferland  and  Goltz,  2001 :46).  This  is 
in  contrast  to  the  flow  pattern  induced  by  the  better-known  groundwater  circulation  wells 
(GCWs)  where  water  circulates  vertically  between  the  two  well  screens  of  a  single  GCW 
(Parsons,  2002). 


Figure  1.  Cross-sectional  view  ofHFTW  operation. 


Direction  of 
Regional 

Groundwater  Flow 


Figure  2.  Plan  view  of  FIFTW  operation. 


HFTWs  offer  several  advantages  over  conventional  treatment  technologies  like 
pump-and-treat  (PAT)  and  penneable  reactive  barriers  (PRBs).  Unlike  PAT  systems, 
FIFTWs  do  not  pump  contaminated  groundwater  to  the  surface.  Eliminating  the  need  to 
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extract  groundwater  is  economically  advantageous  because  it  eliminates  the  expense  of 
(1)  pumping  groundwater  to  the  surface  and  (2)  disposing  of  the  contaminant  (McCarty  et 
al.,  1998:99).  Also,  treating  the  contaminated  water  in  the  subsurface  reduces  risk  of 
exposure  to  the  contaminant  (Gandhi  et  al.,  2002a:4-5).  Another  advantage  of  HFTWs 
when  compared  to  PAT  systems  is  that  after  a  PAT  system  withdraws  water  from  a 
contaminated  zone,  the  hydraulic  gradient  that  was  established  by  the  capture  well 
continues  to  draw  freshwater  into  that  same  zone;  “contaminants  sorbed  to  the  soil  in  that 
zone  desorb  and  contaminate  this  freshwater”  (McCarty  et  al.,  1998:99).  Hence,  pristine 
groundwater  becomes  contaminated  and  unusable.  This  is  especially  problematic  in 
water  shortage  areas,  especially  if  regulations  prohibit  returning  treated  water  to  the 
ground.  HFTWs  are  superior  to  PAT  systems  in  this  regard  because  the  treatment  wells 
inject  water  immediately  after  extracting  it,  so  pristine  groundwater  is  not  drawn  into  the 
contamination  zone  (McCarty  et  al.,  1998:99). 

Recirculation  between  HFTWs  results  in  water  with  dissolved  contaminants 
passing  multiple  times  through  the  treatment  zones,  which  as  noted  earlier  may  either  be 
in-well  or  bioactive  zones  external  to  the  well.  This  recirculation  results  in  higher 
contaminant  destruction  efficiencies  than  are  possible  with  single-pass  systems,  such  as 
PRBs  (Ferland  and  Goltz,  2001:46;  Gandhi  et  al.,  2002:4-5).  Although  PRBs,  like 
HFTWs,  is  an  in  situ  treatment  technology  (Domenico  and  Schwarz,  1998:450),  the  fact 
that  PRBs  are  passive  results  in  significant  disadvantages.  For  one  thing,  PRBs  may  only 
be  applied  under  certain  conditions  (for  instance,  where  the  contamination  plume  is 
relatively  shallow,  so  that  a  barrier  that  can  intercept  the  entire  depth  of  the  plume  can  be 
economically  emplaced).  In  addition,  changing  flow  conditions  can  allow  a 
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contamination  plume  to  bypass  the  PRB  or  result  in  insufficient  residence  time  within  the 
reactive  barrier  to  adequately  treat  the  contaminant.  An  HFTW,  on  the  other  hand,  is  an 
active  technology  that  captures  contaminated  water  by  pumping.  Control  of  flow  and 
recirculation  also  assures  treatment  efficiency  is  adequate. 

The  effectiveness  of  HFTWs  was  demonstrated  at  a  field  test  at  Site  19,  Edwards 
AFB,  CA.  This  field  test  employed  HFTWs  to  remediate  two  aquifers,  an  upper  water 
table  aquifer  and  a  lower  confined  aquifer,  each  contaminated  with  500-1200  pg/L  TCE. 
The  test  involved  pulsed  injections  of  toluene,  oxygen,  and  hydrogen  peroxide  into  a  pair 
of  HFTWs.  The  two  treatment  wells  were  spaced  10  m  apart,  and  each  well  pumped 
contaminated  water  at  a  rate  of  38  L/min  (McCarty  et  al.,  1998:88). 

Each  treatment  well  had  two  screens  (one  in  the  upper  aquifer  and  the  other  in  the 
lower  aquifer)  for  groundwater  to  enter/exit  the  well.  One  treatment  well  extracted 
groundwater  from  the  upper  aquifer,  used  an  in-well  mixer  to  add  toluene,  oxygen,  and 
hydrogen  peroxide  into  the  contaminated  water,  and  discharged  the  water  into  the  lower 
aquifer.  Conversely,  the  second  treatment  well  extracted  groundwater  from  the  lower 
aquifer,  mixed  the  additives  into  the  contaminated  water,  and  discharged  the  water  into 
the  upper  aquifer.  This  strategy  resulted  in  development  of  (1)  in  situ  zones  of 
bioactivity  near  the  discharge  screens  of  the  two  treatment  wells,  where  toluene  was 
aerobically  oxidized,  and  TCE  cometabolized,  and  (2)  groundwater  circulation  cells 
between  the  two  treatment  wells  (McCarty  et  al.,  1998:90).  It  was  demonstrated  that 
HFTW  operation  resulted  in:  (1)  efficient  mixing  of  amendments  into  contaminated 
groundwater  without  the  need  to  pump  the  groundwater  to  the  surface,  and  (2) 
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recirculation  between  the  treatment  wells,  resulting  in  multiple  passes  of  contaminated 
water  through  the  bioactive  zones,  where  treatment  was  effected. 

Regional  groundwater  upgradient  of  the  HFTW  system  had  dissolved  TCE 
concentrations  of  about  1000  jug/L;  downgradient  of  the  system,  groundwater  had  much 
lower  concentrations  of  18  -  24  pg/L,  indicating  97  -  98%  overall  removal  efficiency 
(McCarty  et  al.,  1998:99).  Thanks  to  recirculation  between  the  two  treatment  wells,  it 
was  possible  to  attain  these  high  removals  even  though  only  87%  contaminant  destruction 
was  achieved  with  each  single  pass  of  contaminated  groundwater  through  the  bioactive 
zones  that  were  established  around  the  injection  screens  of  the  treatment  wells  (McCarty 
etal.,  1998:88). 

Christ  et  al.  (1999)  developed  an  analytical  model  to  investigate  application  of 
multiple  injection  and  extraction  well  pairs  to  remediate  groundwater  contaminated  with 
TCE.  Overall  efficiency  (p)  of  an  HFTW  system,  which  is  a  measure  that  compares 
contaminant  concentrations  upgradient  and  downgradient  of  the  system,  is  a  function  of 
the  single-pass  contaminant  removal  efficiency  and  interflow  of  water  between  the 
treatment  wells  (Christ  et  al.,  1999:297,  298).  Christ  et  al.  (1999)  defines  single-pass 
contaminant  removal  efficiency,  psp,  as  “contaminant  removal  for  each  pass  of  water 
through  the  treatment  zone”  (297).  For  the  two-well  HFTW  system  shown  in  Figure  1, 
there  would  be  two  single-pass  contaminant  removal  efficiencies,  pspU  and  pspL,  for  the 
treatment  zones  in  the  upper  and  lower  horizons,  respectively  (Christ,  1997:p3-27). 

To  understand  the  definition  of  interflow,  let  us  assume  we  are  dealing  with  a 
two-well  injection/extraction  system.  In  this  simple  system,  interflow  is  defined  as  the 
water  flowing  into  the  extraction  well  that  originated  in  the  injection  well,  normalized  by 
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the  total  flow  in  the  extraction  well  (Christ  et  al.,  1999:298).  For  a  multi- well  system, 
with  all  wells  pumping  at  the  same  rate,  total  interflow  is  defined  as  the  flow  through  all 
the  extraction  wells  that  originated  in  injection  wells,  normalized  by  the  flow  through  a 
single  extraction  well.  We  may  also  define  a  parameter  that  we  will  call  average 
interflow  (IaVg)  as  the  total  interflow  divided  by  the  number  of  extraction  wells.  Note  for 
the  two-well  system  we  described  above,  total  interflow  equals  average  interflow. 

Interflow  is  a  function  of  system  design  parameters  such  as  pump  rate  and  well 
spacing,  as  well  as  environmental  conditions  such  as  groundwater  regional  flow  and 
aquifer  thickness.  The  reader  may  refer  to  Christ  et  al.  (1999)  for  details  on  how  to 
analytically  calculate  average  interflow  for  given  engineered  and  environmental 
parameters,  under  various  simplifying  conditions  (homogeneity,  horizontal  steady  flow, 
etc.). 

By  definition,  overall  treatment  efficiency  is: 

r|  =  l-Cout/Cin  (2.1) 

where  Qn  and  Cout  are  contaminant  concentrations  upgradient  and  downgradient  of  the 
treatment  system,  respectively  (Christ  et  al.,  1999:304).  By  mass  balance,  overall 
treatment  efficiency  can  also  be  expressed  as 


7  = 


>hP 

l~IaVg(\-Tlsp) 


(2.2) 


If  we  assume  that  the  contaminant  concentration,  Cm,  upgradient  of  a  two-well 
HFTW  system  is  equal  for  both  upper  and  lower  treatment  horizons,  we  may  derive 
expressions  for  concentrations  downgradient  of  the  system  for  both  the  upper  (Coutu)  and 
lower  (C0UtL)  treatment  horizons.  If  we  also  assume  there  is  no  interflow  between  the 
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injection  and  extraction  screens  of  a  single  treatment  well,  the  following  equations  for 


(Coutu)  and  lower  (C0UtL)  apply  (Christ,  1997:p3-27). 


C 


outU 


'(1-^X1-  TJyu  )  + I L  (1  ■ - 1 g  )d  -  rjspL  )(1  -  rjspU  )  N 
v  1  -Vi(1-7^tX1-7Jpu)  , 


(2.3) 


'(1-4  XI  -  VspL )  +  4  (1  ■ -  h  XI  -  yspL  XI  ~  n  spy )' 


(2.4) 


where 


Cin  =  upgradient  concentration,  upper  and  lower  horizons 
Coutu  =  downgradient  concentration,  upper  horizon 
CoutL  =  downgradient  concentration,  lower  horizon 
Iu  =  interflow,  upper  horizon 
II  =  interflow,  lower  horizon 

T)spu  =  single-pass  contaminant  removal  efficiency,  upper  horizon 
fispL  =  single-pass  contaminant  removal  efficiency,  lower  horizon 


Ferland  combined  Christ’s  (1997)  analytical  model  of  groundwater  flow  induced  by 
HFTWs  with  a  submodel  simulating  dehalogenation  of  chlorinated  ethenes  using  a  metal 
catalyst  (Ferland,  2000;  Ferland  and  Goltz,  2001).  The  submodel  simulated  first-order 
destruction  of  chlorinated  ethenes  in  in-well  reactors  containing  palladium-based 
catalysts.  Stoppel  extended  Ferland’s  (2000)  model  to  account  for  catalyst  deactivation 
and  regeneration  (Stoppel,  2001;  Stoppel  and  Goltz,  2002). 

The  disadvantage  to  modeling  HFTWs  using  an  analytical  solution  is  the 
requirement  that  many  simplifying  assumptions  such  as  homogeneity  and  strictly 
horizontal  flow  must  be  made.  Numerical  modeling,  on  the  other  hand,  allows  for  a  more 
accurate  simulation  of  real  conditions,  as  a  numerical  model  can  accommodate 
heterogeneity  and  anisotropy.  Under  these  more  realistic  and  complex  conditions,  a 


22 


numerical  flow  model  like  MODFLOW  (Harbaugh  and  McDonald,  1996)  can  be  used  to 
determine  interflow. 

Both  Huang  and  Goltz  (1998)  and  Gandhi  et  al.,  (2002)  developed  numerical 
flow-and-transport  models  to  simulate  TCE  biodegradation  in  an  HFTW  system.  The  3- 
dimensional  Huang  and  Goltz  (1998)  model  accounts  for  multi-dimensional  flow, 
advective/dispersive  transport  of  dissolved  species,  equilibrium  or  rate-limited  sorption, 
and  biodegradation  (281).  A  partial-implicit  approach  is  used  to  numerically  solve  a  set 
of  nonlinear  partial  differential  equations  describing  fate  and  transport  of  TCE,  toluene, 
oxygen,  and  microbes  (Huang  and  Goltz,  1998:281). 

Gandhi  et  al.  (2002a)  developed  a  3 -dimensional,  numerical  model  to  simulate  the 
HFTW  system  used  at  Site  19,  Edwards  AFB,  which  was  previously  discussed.  Gandhi 
et  al.  (2002a)  used  a  finite-element  approach  in  order  to  1)  enhance  flexibility  in  locating 
both  monitoring  and  treatment  wells  and  to  2)  provide  smaller  grid  dimensions  where 
high  spatial  resolution  is  critical  (i.e.,  near  pumping  wells)  (Gandhi  et  al.,  2002a:  10).  The 
3 -dimensional  model  accounted  for  advective-dispersive  transport,  biodegradation,  and 
inhibition  of  biomass  growth  due  to  hydrogen  peroxide,  which  was  used  as  a  source  of 
oxygen  (Gandhi  et  al.,  2002b:2).  Model  simulations  and  field  data  of  TCE  and  dissolved 
oxygen  concentrations  were  reasonably  consistent  (Gandhi  et  al.,  2002b:26).  Gandhi  et 
al.  (2002b)  points  out  that,  despite  the  heterogeneity  that’s  inherent  in  any  aquifer,  the 
model  adequately  simulated  operation  of  the  HFTW  system  even  when  a  homogeneous 
hydraulic  conductivity  field  was  assumed  (26).  It  appears  that  the  flow  field  imposed  by 
operation  of  the  HFTW  system  results  in  reduction  of  the  effects  of  heterogeneity  on 
treatment  system  perfonnance  (Gandhi  et  al.,  2002b:26). 
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For  this  thesis,  the  system  of  interest  involves  HFTWs  to  remediate  a  perchlorate- 
contaminated  site.  The  author  intends  to  optimize  multiple  objectives  for  this  system,  so 
the  next  section  introduces  multiple-objective  optimization. 


2.4  MULTIPLE  OBJECTIVE  OPTIMIZATION 

A  multi-objective  optimization  problem  (MOP)  consists  of  decision  variables,  two 
or  more  objective  functions,  and  constraints.  These  three  components  of  an  MOP  are 
defined  as  follows: 

•  decision  variables:  variables  whose  numerical  values  are  controlled  by  the 
decision  maker. 

•  objective  function:  a  function  that  maps  decision  variable  values  to  values 
reflecting  a  performance  level;  optimizing  this  performance  level  entails  either 
maximization  or  minimization. 

•  constraints:  restrictions  imposed  by  the  particular  problem  that  must  be  satisfied 
to  render  a  solution  acceptable.  “They  describe  dependences  among  decision 
variables  and  constants  (or  parameters)  involved  in  the  problem”  (Coello  Coello 
et  al.,  2002:5).  If  solutions  satisfy  all  constraints,  then  they  are  feasible; 
otherwise,  they  are  infeasible. 

For  example,  suppose  a  site  engineer  wishes  to  install  a  2-well  pump-and-treat 
(PAT)  system  at  a  site  contaminated  with  TCE.  The  engineer  has  control  over  the 
location  and  pump  rate  of  each  well.  The  engineer  wishes  to  simultaneously  maximize 
contaminant  mass  removal  and  minimize  operating  costs.  The  contaminated  site  is 
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rectangular,  measuring  200  meters  east-west  and  300  meters  north-south.  To  prevent 
excessive  drawdown,  the  pump  rate  of  each  well  must  be  <  200  m3/day. 

In  this  example,  there  are  6  decision  variables:  the  east- west  location  of  each  well 
(ai,  a2),  the  north-south  location  of  each  well  (bi,  b2),  and  the  pump  rate  of  each  well  (Qi, 
Q2).  These  decision  variables  can  be  represented  by  a  single  decision  variable  vector,  x  = 
[ai,  a2,  bi,  b2,  Qi,  Q2].  The  objective  functions  can  be  represented  by  the  vector  f(x)  = 
[fi(x),  f2(x)],  where  fi(x)  =  mass  TCE  removed,  and  f2(x)  =  treatment  cost.  The 
constraints  of  the  problem  are  as  follows: 

0  <  a;  <  200  meters,  for  i  =  1,2 
0  <  bi  <  300  meters,  for  i  =  1,2 
0  <  Qi  <  200  nrVday,  for  i  =  1,2 

Now  consider  a  more  general  formulation  of  a  multiobjective  problem.  Consider 
a  decision  variable  vector  x  =  [xi,  x2. .  .xn]  in  which  elements  xi,  x2. .  ,xn  are  decision 
variables.  Also  consider  an  objective  function  vector  f(x)  =  [fi(x),  f2(x). . .  fk(x)]  in  which 
elements  fi(x),  f2(x). . .  fk(x)  are  objective  functions  to  be  maximized  or  minimized.  Also 
consider  a  constraint  function  vector  g(x)  =  [gi(x),  g2(x)...  gm(x)]  in  which  each  element 
is  a  constraint  of  the  form  g;(x)  <  0  for  i  =  1,  2. .  .m.  Then  an  MOP  is  formally  defined  as 
the  search  for  the  optimum  solution  vector  x*  =  [x*j,  x*2. .  ,x*n]  that  optimizes  f(x)  = 
[fi(x),  f2(x). . .  fk(x)]  and  satisfies  the  constraint  vector  g(x)  =  [gi(x),  g2(x). . .  gm(x)].  The 
constraint  vector  g(x)  establishes  the  feasible  region  Q,  and  any  solution  vector  xefi  is 
defined  as  a  feasible  solution  (Coello  Coello  et  al. ,  2002:7). 
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Unfortunately,  practical  MOPs  do  not  have  a  single  solution  x*eQ  that  optimizes 
all  objective  functions.  Under  nonnal  circumstances,  improving  performance  in  one 
objective  mandates  worsening  performance  in  another;  a  “trade-off’  dilemma  occurs 
among  objectives  that  conflict  with  each  other.  Inability  to  universally  optimize  all 
elements  of  f(x)  renders  the  word  “optimum”  controversial  in  the  MOP  arena.  Therefore, 
determining  the  optimum  solution  ultimately  rests  with  the  preferences  of  the  decision 
maker,  making  the  determination  subjective  rather  than  objective  (Sawaragi  et  al., 
1985:25;  Ringuest,  1992:2). 

Despite  this  subjectivity,  the  concept  of  dominance  enables  the  decision  maker  to 
objectively  distinguish  superior  solutions  from  inferior  ones.  Dominance  is 
mathematically  definable.  Suppose  all  objective  functions  of  f(x)  =  [fi(x),  f>(x). . .  fk(x)] 
are  to  be  minimized  (this  supposition  is  legitimate  because  any  maximization  function 
can  be  expressed  as  a  minimization  function  by  simply  multiplying  the  function  by  -1), 
and  two  vectors  (xi  and  X2  eQ)  exist.  Then,  for  i  =  1,  2. .  .k,  X2  dominates  xi  if  and  only 
if  the  following  condition  is  true  (Coello  Coello  et  al.,  2002: 1 1): 

V  i,  fi(x2)  <  fi(xi)  and  3  i  |  fi(x2)  <  f,(xi)  (2.5) 

The  first  condition  means  that  X2  performs  the  same  or  better  than  xi  with  respect  to  all 
objective  functions;  the  second  condition  means  that  X2  perfonns  better  than  xi  with 
respect  to  at  least  one  objective  function. 

Furthermore,  if  there  is  no  xeQ  that  dominates  X2,  then  X2  is  called 
nondominated,  or  Pareto-optimal.  Still  assuming  all  f(x)  are  to  be  minimized,  the 
vector  x*  is  Pareto-optimal  if,  for  every  xefi  and  i  =  1,  2. .  .k,  either  of  the  following 
conditions  is  true  (Coello  Coello  et  al.,  2002:10): 
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fi(x*)  <  f;  (x)  or 


V  i,  f;(x*)  =  f;  (x) 


(2.6) 


3  i 

Nondomination,  or  Pareto-optimality,  is  an  essential  property  of  any  candidate 
solution  to  an  MOP  (Ringuest,  1992:3).  The  set  of  all  Pareto-optimal  solutions  to  an 
MOP  is  called  the  Pareto-optimal  set,  which  is  defined  as 

P*:=  {xefi  |  — .3  x'eQ,  V  i,  f(x)  <  fi(x')  and  3  i  f,(x)  <  f(x')}  (2.7) 

In  other  words,  each  element  of  P*  is  a  decision  variable  vector  x  that  is  both  feasible  and 
Pareto-optimal.  Coello  Coello  (2002)  emphasizes  that  a  solution’s  membership  in  P* 
depends  on  its  evaluation  with  all  objective  functions  (12-13).  The  mapping  of  all 
solutions  in  P*  to  their  objective  function  values  is  called  the  Pareto  front  (PF*).  The 
Pareto  front  is  defined  as  follows  (Coello  Coello  et  al.,  2002:12): 

PF*:=  {f(x)  =  [fi(x),  f2(x)...  fk(x)]  |  xeP*}  (2.8) 

The  reader  should  understand  that  P*  pertains  to  decision  variable  values  whereas  PF* 
pertains  to  the  objective  function  values.  The  symbolic  relationship  between  P*  and  PF* 
can  be  written  as 

p  $  objective  functions  ^  pp  &  (2  9) 

Each  member  xeP*  is  a  particular  array  of  decision  variable  values.  Each  member 
f(x)ePF*  is  a  particular  combination  of  objective  function  values  corresponding  to  xeP*. 

P*  represents  the  true  Pareto-optimal  set.  In  other  words,  P*  is  an  exhaustive  set 
of  nondominated  solutions  (possibly  infinite).  In  the  real  world,  however,  generating  P* 
is  not  viable  for  several  reasons.  First,  generating  an  infinite  number  of  solutions  is 
impossible,  so  the  solutions  that  are  generated  are  only  a  subset  of  P*.  Second, 
computational  precision  limits  may  render  solutions  that  are  not  identical  to  the  truly 
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Pareto  optimal  solutions.  For  these  reasons,  a  set  of  real-world  MOP  solutions  is  just  a 
discontinuous  approximation  of  P*,  and  this  approximation  is  denoted  Pknown. 

The  distinction  between  “true”  and  “approximate”  also  applies  to  the  Pareto  front. 
For  complex  MOPs,  deriving  an  analytical  expression  for  the  true  Pareto  front  PF*  is 
usually  impossible  (Coello  Coello,  2002:12;  Veldhuizen  and  Lamont,  2000:128);  an 
alternative  method  for  generating  PF*  involves  mapping  each  xeP*  to  its  objective 
function  value  (Coello  Coello,  2002:12).  As  stated  previously,  however,  computational 
limitations  only  permit  us  to  generate  Pknown  instead  of  P*;  so  mapping  Pknown  to  objective 
function  space  yields  PFknown,  which  approximates  PF*.  The  relationship  between 
Pknown  and  PFknown  is  analogous  to  the  relation  ship  between  P*  and  PF*  and  can  be 
represented  symbolically: 


objective  functions 


"  known 


->PF, 


known 


(2.10) 


Given  this  discussion,  the  following  procedure  is  typical  for  solving  real-world 


MOPs: 


•  Generate  Pknown  as  an  estimate  of  P* 

•  Generate  PFknown  by  evaluating  the  objective  function  vector  for  each  xe  Pknown 

•  Choose  a  solution  from  Pknown  that  results  in  an  acceptable  value  of  the 
corresponding  PFknown  on  the  Pareto  front.  “Acceptable”  is  a  subjective  term, 
dependent  upon  the  decision  maker’s  preferences. 

A  popular  stochastic  search  algorithm,  known  as  the  genetic  algorithm  (GA),  is  well- 
adapted  to  solving  complex  MOPs.  More  about  GAs  follows  in  the  next  section. 


28 


2.5  GENETIC  ALGORITHMS 


John  H.  Holland  (1975)  from  the  University  of  Michigan  pioneered  the  stochastic 
search  technique  known  as  the  genetic  algorithm  (GA).  A  genetic  algorithm  searches  for 
solutions  using  principles  analogous  to  evolution  and  biological  genetics.  According  to 
Whitley  (1994),  a  GA  is  “any  population-based  model  that  uses  selection  and 
recombination  operators  to  generate  new  sample  points  in  a  search  space.”  Goldberg 
(1989)  distinguishes  the  genetic  algorithm  from  traditional  search  methods,  like  calculus- 
based,  enumerative,  and  random  search  methods  in  four  ways: 

1 .  GA  requires  that  the  user  encode  the  problem’s  parameters  as  a  data  structure 
(“string”).  Mitchell  (1996)  points  out  that  this  step  is  crucial  for  a  GA  to  be 
successful  (156). 

2.  GA  operates  on  a  population  of  candidate  solutions  -  not  just  one  solution. 

3.  GA  neglects  auxiliary  and  derivative  infonnation;  the  only  information  that  is 
crucial  are  objective  function  values  associated  with  the  strings. 

4.  GA  uses  probabilistic  transition  rules  to  explore  the  search  space  instead  of 
deterministic  rules. 

The  terminology  associated  with  genetic  algorithms  borrows  heavily  from 
biology.  Some  common  tenns  are  defined  here: 

Chromosome:  an  encoded  string  representing  a  complete  parameter  set  or  solution 
to  a  problem;  often  used  interchangeably  with  “string” 

Gene:  a  subunit  of  a  chromosome  dedicated  to  a  particular  parameter 
Allele:  specific  value  that  a  gene  may  assume 
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Genotype:  the  specific  content  of  a  chromosome  that  distinguishes  one 
chromosome  from  another 

Phenotype:  the  decoded  version  of  a  chromosome;  value(s)  obtained  by 
evaluating  a  chromosome  with  its  objective  function(s) 

Fitness:  a  measure  of  how  “good”  a  chromosome  is,  which  is  related  to  its 
evaluation  (phenotype);  a  fitness  function  converts  a  chromosome’s  evaluation 
“into  an  allocation  of  reproductive  opportunities”  (Whitley,  1994), 

Using  the  example  of  the  PAT  system  from  Section  2.4  (multi-objective 
problems),  each  decision  variable  (ai?  a2,  bi,  b2,  Qi,  Q2)  would  be  a  gene ;  the  aggregation 
of  these  variables  is  the  vector  x,  which  a  GA  would  treat  as  a  chromosome.  Using  the 
same  convention  from  section  2.4,  consider  two  hypothetical  vectors  (chromosomes)  xi 
and  X2  such  that 

Xi  =  [156,  187,  275,24,  150,  180]  and  x2  =  [32,  114,  101,299,  100,  200] 
Vectors  Xi  and  x2  consist  of  the  same  decision  variables  {genes),  but  they  assume 
different  values  {alleles).  Because  corresponding  alleles  of  Xi  and  x2  are  dissimilar,  they 
have  different  genotypes.  Suppose  TCE  mass  removal  is  the  only  objective  function  f(x), 
and  f(x0  =  5.2  kg,  and  f(x2)  =  8.7  kg.  The  values  5.2  kg  and  8.7  kg  are  the  phenotypes  of 
xi  and  x2,  respectively.  If  mass  removal  is  the  only  metric  for  chromosome  performance, 
then  chromosome  x2  has  better  fitness  than  xi. 

Goldberg  (1989)  identifies  three  fundamental  GA  operators:  reproduction,  crossover, 
and  mutation  (10).  Reproduction  (or  selection)  means  that  the  GA  reproduces/selects 
chromosomes  from  the  initial  population  and  places  them  into  a  temporary  mating  pool. 
Usually,  the  GA  is  constructed  so  that  the  probability  a  particular  chromosome  is 
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reproduced  is  proportional  to  that  chromosome’s  fitness.  Hence,  chromosomes  with 
higher  fitness  have  higher  probability  of  reproduction;  lower-fitness  chromosomes  have 
less  probability.  This  is  analogous  to  Charles  Darwin’s  notion  of  “survival  of  the  fittest.” 

After  reproduction,  the  mating  pool  is  filled  and  the  chromosomes  are  ready  to 
undergo  crossover.  In  crossover,  the  GA  probabilistically  “mates”  chromosomes  at 
random,  causing  them  to  exchange  portions  of  themselves  with  each  other.  The 
crossover  operator  hopefully  combines  desirable  characteristics  of  those  higher-fitness 
chromosomes  to  make  even  better  chromosomes.  Crossover  probability  usually  ranges 
from  0.4  to  0.9  (Coley,  72:1999).  Finally,  mutation  may  occur.  The  GA,  by  some  small 
probability,  mutates  a  chromosome  by  changing  one  or  more  alleles.  Because  mutation 
probability  tends  to  be  orders  of  magnitude  less  than  crossover  probability  (14),  mutation 
has  less  influence  on  chromosomes  than  crossover.  After  mutation,  the  new  population 
of  chromosomes  is  complete.  Each  cycle  of  fitness  evaluation,  reproduction,  crossover, 
and  mutation  is  called  a  generation. 

The  following  pseudo  code  represents  a  very  simple  GA: 

Initialize  population  of  strings/chromosomes 
Repeat  until  desired  number  of  generations  is  attained 
Evaluate  fitness  (objective  function)  for  each  string 
Reproduce  strings  to  mating  pool  based  on  fitness 
Perform  cross-over 
Mutate  alleles 
End 
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Reproduction,  crossover,  and  mutation  address  the  trade-off  dilemma  between 
exploitation  and  exploration  that  plagues  any  search  algorithm.  Exploitation  is  the 
strategy  of  using  knowledge  of  previous  points  in  the  search  space  to  locate  even  better 
points;  exploration,  on  the  other  hand,  is  the  investigation  of  new  and  unknown  areas  of 
the  search  space  (Beasley  et  al,  1993:63).  Efforts  to  exploit  diminish  efforts  to  explore, 
and  vice  versa.  Although  exploitation  is  useful  because  it  focuses  on  a  promising  area  of 
the  search  space,  excessive  exploitation  causes  the  algorithm  to  neglect  other  areas  of 
search  space  that  may  contain  better  solutions.  Conversely,  exploration  gives  the 
algorithm  freedom  to  “travel”  the  search  space,  but  too  much  exploration  has  the  danger 
of  degenerating  into  a  haphazard  random  search. 

The  reproduction  and  crossover  operators  are  the  exploiting  power  of  the  GA;  they 
work  cooperatively  to  improve  successive  generations  of  chromosomes.  However,  a 
problem  arises  if  these  two  operators  work  strictly  by  themselves:  continuously 
disregarding  lower-fitness  chromosomes  directs  the  search  to  a  very  concentrated  portion 
of  the  search  space  and  diminishes  the  potential  to  find  better  chromosomes  elsewhere. 
This  is  why  the  GA  uses  mutation  to  provide  some  exploring  power.  As  Beasley  et  al 
(1993)  states,  “mutation  provides  a  small  amount  of  random  search,  and  helps  ensure  that 
no  point  in  the  search  space  has  a  zero  probability  of  being  examined”  (60). 

The  reason  GAs  work  is  because  of  the  Fundamental  Theorem  of  Genetic 
Algorithms,  better  known  as  the  Schema  Theorem.  A  schema  (plural  schemata)  is  a 
template  describing  a  subset  of  strings  with  similar  allele  patterns  (Holland,  1992:68; 
Goldberg,  1989:19).  It  is  a  chromosome  template  that  consists  of  alleles  plus  a 
“wildcard”  symbol  (*)  that  can  assume  any  allele. 
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Consider  a  chromosome  population  whose  individual  genes  are  represented  with  the 
binary  alphabet  (0,  1).  The  schema  [11**0*]  has  length  =  6  and  contains  fixed  alleles  in 
the  first,  second,  and  fifth  position;  the  remaining  alleles  are  *’s  that  can  assume  values 
of  0  or  1.  The  schema  [11**0*]  can  represent  the  following  strings:  [110101],  [111000], 
and  [110001]. 

Two  important  characteristics  of  a  particular  schema  is  its  order  and  defining  length. 
Order  is  the  number  of  fixed  alleles  (i.e.  all  alleles  except  *),  and  defining  length  is  the 
distance  between  the  first  and  last  fixed  alleles  (Michalewicz,  1996:46).  The  example 
schema  [11**0*]  has  order  =  3  because  positions  1,  2,  and  5  are  fixed.  Schema  [11**0*] 
has  defining  length  =  4  because  positions  1  and  5  are  the  first  and  last  fixed  alleles  (5  -  1 
=  4). 

The  Schema  Theorem  mathematically  predicts  the  minimum  copies  of  a  specific 
schema  in  the  current  generation  (t)  that  appear  in  the  next  generation  (t  +1).  According 
to  the  Theorem,  schemata  with  1)  short  defining  length  2)  low  order  and  3)  high  fitness 
“receive  exponentially  increasing  trials  in  subsequent  generations”  (Goldberg,  1989:33). 
These  short,  low-order,  high-fitness  schemata  are  called  building  blocks.  An  extension 
of  the  Schema  Theorem  is  the  Building  Block  Hypothesis,  which  states  that  a  GA  “seeks 
near  optimal  perfonnance  through  the  juxtaposition  of  short,  low-order,  high- 
performance  schemata,  or  building  blocks”  (Goldberg,  1989:41).  As  the  GA  repeats  the 
generation  loop,  building  blocks  accumulate  exponenetially  (Holland,  1992:180). 

Although  the  Schema  Theorem  provides  insight,  it  unfortunately  provides  an  inexact 
schema  count  in  the  next  generation  because  it  is  based  on  an  inequality  (>);  hence,  the 
Schema  Theorem  predicts  a  minimum  schema  count  with  no  provision  for  a  maximum, 
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and  it  cannot  accurately  predict  how  a  certain  schema  is  processed  as  generations 
progress  (Whitley,  1994). 

Before  using  a  GA  to  solve  a  problem,  one  must  consider  the  following: 

•  A  GA  does  not  guarantee  a  globally  optimum  solution.  Recall  that  a  GA 
examines  the  search  space  probabilistically,  not  deterministically,  so  the  optimum 
solution  may  be  overlooked. 

•  Constraint  handling  is  a  challenge  for  GAs.  If  a  chromosome  (solution)  violates 
any  constraint,  then  it  is  infeasible  and  theoretically  has  zero  fitness.  However, 
this  situation  causes  problems.  In  practice,  problems  incur  many  constraints,  and 
“finding  a  feasible  point  is  almost  as  difficult  as  finding  the  best”  (Goldberg, 
85:1989).  As  a  result,  the  problem  solver  wants  to  make  use  of  infeasible 
chromosomes  because  they  may  contain  useful  information  (Goldberg,  85:1989; 
Coley,  72:1999).  Goldberg  (1989)  identifies  the  “penalty  method”  as  a  means  of 
retaining  infeasible  solutions.  In  the  penalty  method,  an  infeasible  solution  is 
allowed  to  have  fitness,  but  its  fitness  degrades  proportionally  to  its  constraint 
violation  (85).  As  Coley  (1999)  points  out,  “the  form  of  the  penalty  function  must 
be  chosen  with  care”  (72). 

•  According  to  Beasley  (1993),  “ convergence  is  the  progression  towards  increasing 
uniformity”  (60),  and  premature  convergence  is  a  classic  problem  with  GAs.  In 
premature  convergence,  individuals  with  relatively  high,  but  not  the  best,  fitness 
rapidly  dominate  the  population,  causing  the  population  to  converge  to  a  local 
maximum.  As  virtually  identical  strings  continue  reproduction  and  crossover,  the 
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population  becomes  homogeneous,  and  continued  search  efforts  become  futile 
(Beasley,  64,  1993). 

•  Fitness  evaluation  of  all  chromosomes  is  undoubtedly  the  most  time-intensive 
portion  of  the  GA.  Because  a  GA  operates  on  a  population  of  candidate  solutions, 
it  incurs  the  burden  of  evaluating  the  entire  population;  therefore,  the  evaluation 
function  “must  also  be  relatively  fast”  (Whitley,  1994).  Extreme  computational 
effort  may  necessitate  networking  processor  hardware  (Coley,  84-85:1999). 
Assigning  sensible  values  for  internal  GA  parameters  like  crossover  probability,  mutation 
probability,  and  generation  size  is  a  dubious  task  because  there  are  “no  conclusive  results 
on  what  is  best”  (Mitchell,  1996:175).  Appropriate  parameter  values  depend  upon  the 
nature  of  the  problem  (Coley,  22:1999). 


2.6  PARETO-BASED  MULTI-OBJECTIVE  GENETIC  ALGORITHMS 

This  section  ties  together  sections  2.4  and  2.5  to  introduce  the  multi-objective 
genetic  algorithm  (MOGA).  In  particular,  this  section  focuses  on  MOGAs  that  use 
Pareto-based  approaches. 

The  reader  is  reminded  that  all  MOGAs  are  stochastic  search  algorithms  based 
upon  probabilistic  methods.  Because  they  are  not  deterministic,  MOGAs  can  only 
generate  an  approximated,  incomplete  version  of  P*  rather  than  P*  itself.  Therefore,  the 
goal  of  a  Pareto-based  MOGA  is  convergence  of  Pknown  towards  P*.  MOGAs,  like 
single-objective  GAs,  operate  on  a  population  of  candidate  solutions  (chromosomes)  as 
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opposed  to  a  single  solution;  therefore,  the  strength  of  an  MOGA  is  its  ability  to  uncover 
multiple  nondominated  solutions  (Pknown)  in  a  single  run. 

MOGAs  differ  from  single-objective  GAs  in  two  fundamental  ways: 

•  Solutions  quantity.  A  single-objective  GA  typically  searches  for  a  single  solution 
that  best  optimizes  a  single  objective  function.  A  MOGA,  however,  searches  for 
multiple  solutions  that  are  Pareto-optimal. 

•  Fitness  assignment.  A  single-objective  GA  assigns  higher  fitness  to  solutions  that 
yield  better  performance  relative  to  a  single  objective.  However,  MOGAs  operate 
on  two  or  more  conflicting  objectives;  therefore,  Pareto-dominance  is  the  only 
basis  for  assigning  fitness  (Fonseca  and  Fleming,  1995:46). 

2.6.1  Sorting  nondominated  solutions 

Solutions  exhibit  various  degrees  of  nondomination  relative  to  each  other.  For 
example,  for  any  two  feasible  solutions  xi  and  x2,  the  following  scenarios  are  possible: 

•  Solution  xi  dominates  x2 

•  Solution  x2  dominates  xi 

•  Neither  solution  dominates  the  other 

Because  nondomination  typically  influences  selection  in  MOGAs,  there  is  a  need  to  sort, 
or  rank,  chromosomes  according  to  their  degree  of  nondomination.  Veldhuizen  and 
Lamont  (2000)  point  out  that  MOGAs  typically  rely  on  the  ranking  schemes  of  Goldberg 
(1989)  and  Fonseca  and  Fleming  (1993).  With  the  Goldberg  (1989)  method,  the  MOGA 
initially  identifies  nondominated  chromosomes  within  the  current  population,  assigns 
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them  rank  =  1,  and  removes  them  from  the  population.  Next,  the  MOGA  identifies 
nondominated  chromosomes  in  the  remaining  population,  assigns  them  rank  =  2,  and 
removes  them  from  the  population.  This  cycle  of  identification,  ranking,  and  removal 
continues  until  all  chromosomes  in  the  current  population  are  ranked.  Reproductive 
probabilities  are  then  based  on  rank  (Goldberg,  1989:201). 

The  Fonseca  and  Fleming  (1993)  method  is  somewhat  different.  Suppose  x  is  a 
particular  chromosome  in  the  current  population.  Then  x  receives  rank  =  1+  p,  where  p  is 
the  number  of  chromosomes  in  the  entire  current  population  that  dominate  x.  Like 
Goldberg’s  (1989)  method,  the  Fonseca  and  Fleming  (1993)  method  ensures  all 
nondominated  chromosomes  in  the  current  population  have  rank  =  1 .  However,  while 
Goldberg’s  (1989)  method  successively  removes  individuals  after  ranking  them,  the 
Fonseca  and  Fleming  (1993)  method  retains  all  chromosomes  in  the  current  population  to 
be  ranked  among  each  other  simultaneously. 

2.6.2  Maintaining  diverse  solutions 

In  solving  MOPs,  Pareto-optimality  is  an  essential  property  of  any  solution,  but 
this  alone  is  not  sufficient.  A  secondary  requirement  is  that  solutions  reflect  sufficient 
diversity.  Diversity  can  apply  to  either  the  solutions  themselves  (P*)  or  to  their 
evaluation  (PF*).  This  discussion,  however,  assumes  that  solution  phenotype  is  more 
important  than  solution  genotype,  and,  therefore,  limits  discussion  of  diversity  to  PF*, 
which  reflects  how  solutions  perform  with  the  objective  functions.  A  MOGA  may 
generate  thousands,  or  even  millions,  of  Pareto-optimal  solutions,  but  if  they  all  converge 
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to  the  same  portion  of  PF*,  the  large  solution  quantity  yields  no  benefit  for  the  decision 
maker. 

Pareto-based  MOGAs  assign  nondominated  chromosomes  equal  fitness,  which 
theoretically  means  they  all  have  equal  opportunity  for  reproduction.  In  practice, 
however,  this  is  not  the  case  because  a  GA  operates  on  a  finite  population  as  opposed  to 
an  infinite  population.  Repetitious  sampling  of  small,  finite  populations  causes  stochastic 
errors  to  accumulate.  Consequently,  even  when  multiple  solutions  offer  no  advantage 
over  each  other,  the  population  ultimately  converges  to  a  single  solution  (Goldberg  & 
Richardson,  1987:41-42).  Such  convergence  due  to  stochastic  errors  from  repeatedly 
sampling  small  population  sizes  is  called  genetic  drift  (Goldberg,  1989:185-186). 

Genetic  drift  causes  MOGAs  to  eventually  converge  to  a  single  nondominated  solution 
(Fonseca  and  Fleming,  1993:6/2),  which  results  in  “crowding”  (convergence)  on  the 
Pareto-front  and  loss  of  diversity.  Goldberg  and  Richardson  (1987)  demonstrated  genetic 
drift  with  single-objective  GA  optimizing  a  function  with  multiple  optima  or  “peaks” 

(i.e.,  multimodal  function).  To  counteract  genetic  drift,  Goldberg  and  Richardson  (1987) 
developed  the  idea  of  fitness  sharing. 

Fitness  sharing  is  fundamentally  driven  by  a  parameter  called  the  niche  radius, 
share-  The  parameter  ashare  is  a  user-defined,  radial  distance  that  defines  a  circular  niche 
around  each  point  in  objective  space.  Points  that  lie  within  ashare  of  each  other  “share” 
the  same  niche  and  degrade  each  other’s  fitness  (i.e.  they  must  “share”  each  other’s 
fitness).  To  determine  the  degraded/shared  fitness  of  chromosome  X;  due  to  nearby 
points  in  its  niche,  a  sharing  function  is  used.  Goldberg  and  Richardson  (1987)  identify 
the  following  sharing  function: 
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Sh(dij)  =  1  -  (dij  / CJ share)01  if  dij  <  Gshare  (2.11) 

=  0  otherwise 

Equation  2.1 1  depends  upon  parameters  ashare  and  a.  The  argument  dij  is  the  distance  in 
objective  space  between  \,  and  some  other  chromosome  xj;  Xj  can  be  a  member  of  the 
entire  population,  a  member  of  the  same  equivalence  class  as  X;,  or  some  other  user- 
defined  subset.  The  niche  count  m,  of  chromosome  x;  is  the  summation  of  its  sharing 
function  values  with  each  chromosome  xj: 

ffl,  =  E*  SMd,,)  (2.12) 

The  niche  count  m ,  quantifies  the  magnitude  of  crowding  around  chromosome  x;  in 
objective  space.  Finally,  the  fitness  of  Xj  is  degraded  (“shared”)  by  dividing  its  original 
fitness  f  by  its  niche  count  m,;  shared  fitness  =  f  /  m,.  Probabilities  for  selection  are  then 
based  upon  the  shared  fitness  values  (Goldberg  and  Richardson,  1987). 

In  MOGAs,  fitness  sharing  can  resolve  the  dilemma  of  how  to  select  equal¬ 
ranking  chromosomes.  For  example,  a  MOGA  may  initially  assign  equal  fitness  values 
to  all  chromosomes  of  some  rank  =  k.  To  encourage  diversity,  the  MOGA  reduces  fitness 
of  each  Gran  king  chromosome  in  proportion  to  its  niche  count.  These  adjusted  fitness 
values  then  guide  selection.  It  is  also  possible  to  implement  fitness  sharing  without 
directly  adjusting  chromosome  fitness  (Horn  et  a/.,  1994). 

2.6.3  Mating  restriction. 

Arbitrary  recombination  of  chromosome  pairs  having  extreme  dissimilarities 
(genotypic  or  phenotypic)  may  result  in  low-performance  offspring  called  “lethals” 
(Goldberg,  1989:184).  Therefore,  it  makes  sense  to  restrict  crossover  to  chromosome 
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pairs  that  are  similar.  As  Goldberg  (1989)  points  out,  the  sensibility  of  this  restriction 
follows  from  nature,  which  prohibits  mating  among  different  organisms  (mammals,  birds, 
reptiles,  etc.).  Restriction  of  crossover  among  chromosome  pairs  is  called  mating 
restriction. 

The  purpose  of  mating  restriction  in  MOGAs  is  to  prevent  creation  of 
chromosomes  that  are  not  Pareto-optimal.  Mating  restriction  within  MOGAs  is  more 
relevant  to  phenotype  than  genotype  due  to  greater  interest  in  generating  the  Pareto  front 
(Veldhuizen  and  Lamont,  2000).  The  distance  parameter  amate  controls  mating 
restriction;  chromosome  x;  cannot  crossover  with  other  chromosomes  that  reside  more 
than  omate  from  x,  in  objective  space.  If  a  MOGA  uses  both  fitness  sharing  and  mating 
restriction,  the  common  practice  is  to  set  oshare  =  omate,  although  such  practice  is  purely 
arbitrary.  No  sound  theory  justifies  the  inclusion/exclusion  of  mating  restriction  in 
MOGAs  (Veldhuizen  and  Lamont,  2000). 

2.6.4  Pareto-based  MOGAs  in  the  literature. 

Horn  et  al.  (1994)  developed  the  niched  Pareto  genetic  algorithm  (NPGA),  which 
uses  Pareto  domination  tournaments  and  fitness  sharing  chromosome  selection.  The 
NPGA  creates  a  tournament  by  randomly  picking  two  solutions  Xi  and  x2  from  the 
current  population  and  by  randomly  filling  a  secondary  population  (the  comparison  set) 
with  solutions;  the  comparison  set  is  so-named  because  both  xi  and  x2  are  compared  to 
each  member  of  the  comparison  set.  If  the  comparison  set  dominates  xi  but  not  x2,  then 
x2  “wins  the  tournament”  and  is  selected  for  reproduction.  If  the  comparison  set  either 
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(a)  dominates  both  xi  and  x2  or  (b)  dominates  neither  xi  nor  X2,  then  the  tournament  ends 
in  a  tie.  In  the  case  of  a  tie,  it  is  probable  that  xi  and  x2  belong  to  the  same  equivalence 
class,  or  partial  order.  The  NPGA  resolves  ties  with  equivalence  class  sharing.  In 
equivalence  class  sharing,  the  NPGA  plots  each  solution  in  objective  space  and  uses  the 
user-defined  parameter  Gshare  to  assess  crowding  of  adjacent  points.  Horn  et  al.  (1994), 
however,  defines  niche  count  differently  than  equation  2.12.  With  the  NPGA,  the  niche 
count  nii  of  chromosome  x;  simply  equals  the  count  of  points  within  ashare  of  x;  in 
objective  space.  The  NPGA  ultimately  selects  either  xi  or  x2  based  on  lowest  niche  count 
(i.e.  least  crowded  in  objective  space). 

Srinivas  and  Deb  (1994)  introduced  the  nondominated  sorting  genetic  algorithm 
(NSGA),  so-named  because  it  is  based  on  a  nondominated  sorting  procedure.  The  key 
features  of  the  NSGA  are  1)  a  ranking  selection  method  to  emphasize  favorable  solutions 
and  2)  fitness  sharing  to  form  niches.  The  NSGA  first  identifies  nondominated 
chromosomes  in  the  current  population  and  assigns  them  rank  =  1 .  The  NSGA  initially 
assigns  the  same  “dummy”  fitness  value  to  these  rank-1  chromosomes.  The  NSGA 
employs  fitness  sharing  by  dividing  each  rank-1  chromosome’s  fitness  by  its  niche  count, 
calculated  with  equations  2.11  (a=2)  and  2.12.  The  NSGA  then  separates  the  rank- 1 
chromosomes  from  the  population,  identifies  nondominated  chromosomes  in  the 
remaining  population,  assigns  them  the  next  rank  (rank  =  2),  and  gives  them  all  the  same 
fitness  value,  which  is  less  than  the  minimum  shared  fitness  of  the  previous  rank  (or 
“front”).  Hence,  the  NSGA  uses  the  Goldberg  (1989)  ranking  method.  The  process  of 
identifying  nondominated  solutions,  assigning  next  rank,  and  adjusting  fitness  repeats  for 
all  chromosomes.  The  adjusted  fitness  values  determine  selection  probabilities. 


41 


Deb  et  al.  (2002)  developed  the  Non-dominated  Sorting  GA-II  (NSGA-II)  to 
rectify  flaws  of  the  original  NSGA  of  Srinivas  and  Deb  (1994).  Deb  et  al.  (2002) 
identifies  three  key  flaws  with  NSGA,  which  are  1)  high  computational  complexity  of 
nondominated  sorting,  2)  lack  of  elitism,  and  3)  reliance  on  user-specified  Gshare 

To  improve  nondominated  sorting,  NSGA-II  compares  each  solution  in  the 
current  population  P  with  a  partially-filled  population  (P')  of  nondominated  solutions 
rather  than  the  entire  P.  Comparing  members  of  P  with  P'  instead  of  the  entire  P  itself 
has  good  potential  for  reducing  the  number  of  comparisons  that  the  algorithm  must 
check.  NSGA-II,  like  NSGA,  uses  the  Goldberg  (1989)  ranking  method  with 
nondominated  fronts  sequentially  ranked  and  removed  from  P. 

Deb  et  al.  (2002)  criticized  NSGA’s  fitness  sharing  because  of  its  reliance  on  the 
user  chosen  parameter  Gshare-  Computational  burden  is  apparent  in  equations  2.11  and 
2.12;  NSGA  must  calculate  phenotypic  distances  between  chromosome  x,  and  every  Xj  in 
the  same  front  and  then  evaluate  the  sharing  function  for  all  these  distances.  Also,  the 
utility  of  equation  2.1 1  is  sensitive  to  Gshare  chosen  by  the  user. 

To  overcome  these  problems,  Deb  et  al.  (2002)  developed  a  fitness-sharing 
procedure  that  minimized  the  need  for  distance  calculations  and  did  not  require 
specification  of  a  niche  radius  parameter.  NSGA-II  sorts  the  population  in  ascending 
order  using  each  series  of  objective  function  values.  Boundary  points  (minimum  and 
maximum)  are  assigned  infinite  crowding  distance  since  they  lack  an  adjacent  neighbor. 
For  each  intermediate  point  i,  NSGA-II  computes  only  two  crowding  distances:  the 
distances  between  i  and  its  two  adjacent  neighbors.  The  overall  crowding  distance  value 
of  point  i  equals  the  sum  of  individual  distance  values  associated  with  each  objective. 
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The  disadvantage  of  this  procedure  is  the  computational  burden  of  sorting  points  for  each 
series  of  objective  function  values  (Deb  et  al.,  2002). 

The  main  loop  of  NSGA-II  works  as  follows:  NSGA-II  combines  the  parent  and 
children  populations  and  sorts/ranks  the  combined  population  using  the  improved 
nondominated  sorting  procedure  explained  earlier.  Next,  the  algorithm  computes 
crowding  distances  for  each  front.  NSGA-II  fdls  an  A-size  mating  pool  with  priority 
given  to  lower-ranking  chromosomes;  in  this  regard,  NSGA-II  uses  an  elitist  method. 
When  the  algorithm  encounters  a  particular  front/rank  whose  members  exceed  space  left 
in  the  pool,  the  algorithm  resorts  to  fitness  sharing.  NSGA-II  picks  chromosomes  from 
the  front  that  have  the  lowest  crowding  distances.  The  filled  mating  pool  undergoes 
crossover  and  mutation  to  create  a  new  population  of  children  of  size  A  (Deb  et  al., 
2002). 

Zitzler  and  Thiele  (1998)  developed  the  strength  Pareto  evolutionary  algorithm 
(SPEA)  which  combines  traditional  and  new  MOGA  techniques  to  find  nondominated 
solutions.  SPEA  starts  by  initializing  a  population  ( P )  and  an  external  dominated  set 
( P  ),  which  is  initially  empty.  The  algorithm  copies  nondominated  chromosomes  from  P 
to  P  'and  removes  chromosomes  from  P  'that  are  dominated  by  any  other  member  of  P '. 
A  specified  size  parameter  (A)  limits  the  number  of  nondominated  chromosomes  that  P' 
can  store,  and  a  clustering  technique  known  as  the  average  linkage  method  is  used  to 
“prune”  (remove)  excess  chromosomes  from  P '.  SPEA  then  computes  “strength” 
(fitness)  for  each  x;  e  P '  by  computing  a  strength  value  (s,)  that  is  directly  proportional 
to  the  number  of  chromosomes  in  P  'that  x,  dominates.  SPEA  also  computes  a  strength 
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value  (.s',)  for  each  x,  e  P  by  summing  the  strength  values  of  all  chromosomes  in  P  'that 
dominate  Xj.  Zitzler  and  Thiele  (1998)  reverse  the  meaning  of  “fitness”  within  SPEA 
because  higher  reproductive  probabilities  correspond  to  lower  fitness/strength  as  opposed 
to  higher  fitness.  Furthermore,  the  SPEA  niching  technique  defines  niches  in  tenns  of 
Pareto  dominance  instead  of  phenotypic  distance,  which  eliminates  the  need  for  a  user- 
specified  distance  parameter  (i.e.  ashare).  SPEA  uses  binary  tournament  selection  with 
replacement  to  copy  chromosomes  from  both  populations  ( P  u  E)  to  the  mating  pool. 
Once  full,  chromosomes  in  the  mating  pool  undergo  problem-specific  crossover  and 
mutation,  which  completes  a  generation. 

Van  Veldhuizen  (1999)  developed  the  “messy”  multiobjective  genetic  algorithm 
(MOMGA),  which  differs  from  the  previous  MOGAs  in  that  MOMGA  principally 
operates  on  building  blocks  (BBs)  (see  Section  2.5)  instead  of  whole  chromosomes. 
MOMGA  consists  of  three  distinct  phases:  (1)  initialization,  (2)  primordial  phase,  and  (3) 
juxtapositional  phase.  In  the  initialization  phase,  MOMGA  uses  Partially  Enumerative 
Initialization  (PEI)  to  initialize  a  population  representing  all  possible  BB  variations  of 
specific  size.  Following  PEI,  building  blocks  are  evaluated  with  respect  to  a  competitive 
template,  and  then  MOMGA  performs  subroutines  for  the  primordial  and  juxtapositional 
phases.  In  the  primordial  phase  subroutine,  Pareto-based  tournament  selection  repeatedly 
adjusts  the  size  of  the  current  BB  population.  In  the  juxtapositional  phase,  MOMGA 
clones  desired  BBs  and  recombines  them  into  complete  chromosomes  using  a  “cut-and- 
splice”  operator.  New  chromosomes  are  subsequently  evaluated  against  the  competitive 
template  and,  if  they  are  currently  Pareto  optimal,  they  are  added  to  the  currently  Pareto 
optimal  set.  Once  the  juxtapositional  loop  ends,  MOMGA  updates  the  competitive 
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template  using  the  best  fitness  values.  MOMGA  then  executes  another  cycle  of 
initialization,  primordial,  and  juxtapositional  phases  (Van  Veldhuizen,  1999). 

The  disadvantage  of  MOMGA  lies  in  the  PEI,  which  deterministically  enumerates 
all  possible  BB  clones  and,  consequently,  incurs  significant  computational  expense 
(Zydallis  et  al.,  2001).  MOMGA-II,  developed  by  Zydallis  et  al.  (2001),  is  similar  to 
MOMGA  with  some  notable  exceptions.  MOMGA-II  uses  Probabilistically  Complete 
Initialization  (PCI)  instead  of  PEI  to  initialize  a  limited-size  BB  population,  which 
effectively  reduces  computational  “bottlenecks.”  Also,  instead  of  a  primordial  phase, 
MOMGA-II  implements  a  Building  Block  Filtering  (BBF)  phase  to  probabilistically 
ensure  that  all  desirable  BBs  are  in  the  initial  population.  BBF  essentially  reduces  the 
number  of  BBs  and  stores  the  best  BBs  found  (Zydallis  et  al.,  2001). 
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3.0  METHODOLOGY 


3.1  INTRODUCTION 

In  this  chapter,  we  describe  the  procedures  used  to  achieve  our  research  objective 
of  determining  values  of  design  variables  that  optimize  technology  cost  and  performance 
under  various  contaminated-site  conditions.  In  Section  3.2  of  this  chapter,  we  present  the 
details  of  Parr’s  (2002)  technology  model  for  remediating  perchlorate  contamination 
using  HFTWs  in  conjunction  with  in  situ  biodegradation.  In  Section  3.3,  we  establish 
parameter  values  to  describe  an  actual  perchlorate-contaminated  site.  In  Section  3.4,  we 
formulate  a  multi-objective  optimization  problem,  establishing  objective  functions  for 
application  of  the  remediation  technology  at  the  site.  In  Section  3.5  we  present  a  multi¬ 
objective  GA  (MOGA)  that  is  used  to  solve  the  optimization  problem.  Finally,  in  Section 
3.6,  we  describe  the  details  of  implementing  the  MOGA. 


3.2  TECHNOLOGY  MODEL 

Parr  (2002)  developed  a  technology  model  that  combines  groundwater  flow  and 
transport  of  dissolved  species  induced  by  operation  of  an  HFTW  system  with  a 
bioremediation  submodel  that  simulates  perchlorate  reduction  due  to  introduction  of  an 
electron  donor  in  the  HFTW  treatment  wells.  In  this  section,  we  describe  both  the  flow- 
and-transport  model  and  the  bioremediation  submodel  implemented  in  this  study. 
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3.2.1  Flow  and  Transport  Model 


Parr  (2002)  adapted  the  Huang  and  Goltz  (1998)  numerical  model,  which  was  originally 
used  for  TCE,  to  simulate  flow  and  transport  of  CIO4"  in  groundwater  flow  fields  induced 
by  a  2-well  HFTW  system.  The  3-D  flow-and-transport  model  is  a  numerical  model  that 
simulates  advection,  dispersion,  and  consumption  of  four  individual,  dissolved  species 
(electron  donor,  oxygen,  nitrate,  and  perchlorate).  Consumption  rates  of  dissolved 
species  are  due  to  microbially  mediated  redox  reactions,  and  rate  equations  for  these 
reactions  are  in  the  biological  submodel  (section  3.2.2).  The  flow-and-transport  model 
assumes  only  the  electron  donor  (acetate  CH3COO')  can  sorb,  and  this  sorption  is 
assumed  to  be  an  equilibrium  process  that  is  both  linear  and  reversible.  The  remaining 
species  (O2,  NO3",  CIO4’)  are  non-sorbing.  Parr  (2002)  further  assumes  microorganisms 
remain  stationary  (i.e.  fixed  to  the  aquifer  material).  The  flow-and-transport  model 
equations  are  presented  in  Appendix  A. 

The  flow-and-transport  model  is  coded  in  FORTRAN.  MODFFOW  (Harbaugh 
and  McDonald,  1996),  which  is  a  three-dimensional  finite  difference  code,  is  used  as  a 
subroutine  to  calculate  groundwater  flow,  as  the  model  transport  equations  (A.1-A.4) 
require  groundwater  flow  velocities  (v)  throughout  the  problem  domain  as  input.  Given 
boundary  conditions  for  hydraulic  head  (h),  along  with  values  for  aquifer  hydraulic 
conductivity  (K),  porosity  (n),  HFTW  treatment  well  locations,  and  treatment  well 
pumping  rates,  MODFFOW  solves  the  steady-state  flow  equation  (3.1)  to  determine  the 
hydraulic  head  field  and  then  applies  Darcy’s  law  (3.2)  to  compute  the  steady-state 
velocity  field  in  three  dimensions. 
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W2h  =  0 


(3.1) 


v  =  — ——  V  h  (3.2) 

n 

The  flow-and-transport  model  separately  accepts  other  input  information,  like 
injected  donor  concentration  and  injection  frequency.  After  accepting  the  three 
dimensional  steady-state  groundwater  flow  velocity  input  from  MODFLOW,  along  with 
initial/boundary  concentrations  of  the  dissolved  species  and  microbes,  the  flow-and- 
transport  model  uses  a  self-adaptive,  partial  implicit  finite  difference  technique  to 
numerically  solve  the  partial  differential  equations  describing  advection,  dispersion, 
sorption,  and  consumption  in  Appendix  A.  The  model  can  provide  concentrations  for  all 
species  (electron  donor,  electron  acceptors,  and  microorganisms)  at  any  location  over 
time. 

3.2.2  Biological  Treatment  Submodel 

Parr  (2002)  based  his  biological  treatment  submodel  on  a  perchlorate 
biodegradation  model  developed  by  the  environmental  firm  Envirogen.  The  biological 
treatment  submodel  model  consists  of  differential  equations,  based  on  dual-Monod 
kinetics,  that  describe  1)  consumption  rates  of  dissolved  species  (electron  donor,  oxygen, 
nitrate,  perchlorate)  due  to  microbial  redox  reactions,  and  2)  biomass  changes.  A  notable 
feature  of  the  submodel  is  that  it  addresses  the  competition  among  multiple  electron 
acceptors  for  oxidation  of  the  electron  donor;  the  submodel  assumes  that  oxygen  is 
preferentially  reduced  over  nitrate,  which  is  preferentially  reduced  over  perchlorate.  The 
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reader  may  refer  to  Appendix  A  for  a  more  detailed  discussion  of  the  submodel 
equations. 


3.3  SITE  MODEL 

Applying  a  technology  model  to  simulate  remediation  at  an  actual  contaminated 
site  necessitates  translating  the  site  characteristics  into  model  parameters.  Using  a  pre¬ 
processor  such  as  the  one  available  in  Visual  MODFLOW,  a  grid  composed  of  discrete, 
rectangular  cells  can  be  easily  constructed  to  approximate  a  perchlorate-contaminated 
site.  A  key  advantage  of  the  Visual  MODFLOW  pre -processor  and  the  MODFLOW 
code  is  that  they  allow  the  user  to  represent  the  site  with  non-uniform  cell  sizes.  The  user 
can  specify  smaller,  numerous  cells  for  better  resolution  in  the  vicinity  of  pumping  wells 
and  treatment  wells;  and  the  user  can  conversely  specify  larger,  fewer  cells  in  less  critical 
areas  to  reduce  computational  time.  Using  the  pre-processor,  the  user  can  locate 
pumping  and  injection  wells  anywhere  in  the  finite  difference  grid,  establish  head  and 
concentration  boundary  conditions,  contaminant  concentration  initial  conditions,  and 
specify  hydrogeological  characteristics.  Site  information  is  saved  to  data  files  that  are 
later  read  by  MODFLOW. 

A  plan  view  of  the  site  is  shown  in  Figure  3.1.  This  finite  difference,  composite 
grid  is  actually  3 -dimensional,  185  m  long  by  185  m  wide  by  32  m  deep.  This  grid  is  a 
composite  of  discrete,  rectangular  cells  of  non-uniform  size.  Individual  cells  have  sides 
of  3  m,  5  m,  or  10  m.  More  numerous,  smaller  cells  (3  m  x  3  m)  are  located  within  the 
vicinity  of  the  pumping  wells  for  better  resolution.  Larger  cells  (10  m  x  10  m)  are 


49 


specified  at  the  extremities  of  the  grid,  where  resolution  is  less  important,  to  reduce 
computational  time.  This  grid  has  four  layers  with  a  unifonn  horizontal  hydraulic 
conductivity  that  is  twenty  times  greater  than  the  vertical  conductivity.  This  anisotropy  is 
assumed  constant  over  the  entire  modeled  volume.  The  top  layer  represents  an  8  meter 
deep  zone,  where  the  water  table  is  located  an  average  of  1.5  meters  below  the  surface. 
The  second  and  fourth  layers  (10  meters  deep  each)  are  where  the  upper  and  lower 
screens  of  the  treatment  wells  are  located,  and  the  third  layer  (4  meters  deep)  separates 
the  screened  intervals. 

Two  treatment  wells  are  oriented  perpendicular  to  the  direction  of  groundwater 
flow,  which  is  left  to  right  in  Figure  3.1.  The  configuration  in  Figure  3. 1  is  used  to 
evaluate  how  well  a  pair  of  HFTWs  can  be  used  to  contain  a  perchlorate  plume  from 
being  transported  downgradient  from  a  continuous  source. 

MONITORING 

WELLS 

I 


Figure  3.  Contaminated  site  model  -  aerial  view. 
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3.4  FORMULATION  OF  MULTI-OBJECTIVE  PROBLEM 


Having  developed  a  technology  model  and  a  site  model,  we  are  now  in  a  position 
to  fonnulate  a  multi-objective  problem.  Over  a  given  time  of  technology  operation,  we 
wish  to  1)  destroy  as  much  perchlorate  as  possible  and  simultaneously  2)  keep  operating 
costs  of  the  remediation  technology  low.  These  objectives  can  be  represented  as  fi  and 
fj,  respectively.  Next,  we  must  identify  characteristics  (decision  variables)  that  a  site 
engineer  could  manipulate  to  pursue  these  objectives.  Such  characteristics  are 
Q  =  pump  rate  (m  /day)  for  each  well  in  the  HFTW  well  pair 
d  =  spacing  between  the  two  treatment  wells  in  the  well  pair  (meters) 

Cm  =  injected  concentration  of  acetate  (mg/L) 
p  =  acetate  injection  pulse  duration  (in  32nds  of  a  day) 

We  now  want  to  mathematically  express  the  objectives  fi  and  fz  as  functions  of  the 
decision  variables.  However,  objective  fi  (mass  perchlorate  destroyed)  cannot  be 
explicitly  written  as  a  function  of  the  decision  variables.  Looking  at  the  technology 
model,  we  can  see  that  fi  is  a  function  of  the  decision  variables,  but  determining  fi  for  a 
given  set  of  decision  variable  values  requires  numerical  evaluation  of  the  set  of  partial 
differential  equations  that  comprise  the  technology  model.  Therefore,  we  generically 
represent  objective  fi  as 

Mass  CKV  destroyed  =  U{Q,  d,  Cin,p ) 

On  the  other  hand,  the  relationship  between  f2  (total  operating  cost)  and  [ Q ,  d,  Cin,  p]  can 
be  explicitly  formulated.  To  simplify  our  comparison  of  operating  costs  for  different 
implementations  of  the  remediation  technology,  we  assume  that  operating  cost 
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differences  are  only  due  to  differences  in  the  a)  cost  of  electron  donor  and  b)  cost  of 
operating  the  pumps.  That  is,  we  implicitly  assume  that  capital  costs,  as  well  as  other 
recurring  costs  (e.g.  maintenance)  for  different  technology  implementations  are  equal. 

Cost  of  electron  donor  depends  upon  how  much  electron  donor  is  injected  over 
the  duration  of  the  remediation  period  t,  and  we  can  write  this  cost  as 

Material  cost  =  2  *  Q*  t*  C,n  *  p/32  *  1000  L/m3  *  Price  donor,  where 
2  =  number  of  treatment  wells 
t  =  treatment  period  (days) 

32  =  maximum  pulse  duration  (32  pulse  units  =  1  day) 

1000  L/m'  =  conversion  factor 

Price  donor  =  price  of  electron  donor  injected  ($/mg  donor) 

The  remediation  technology  also  incurs  the  cost  of  operating  the  pumps  in  the  two 
HFTWs.  Assuming  continuous  pump  operation,  the  pump  cost  equation  becomes 

Pump  cost  =  2  *  Q  *  t  *  E  *  Price eiec 

where 

2  =  number  of  treatment  wells 

-3 

E  =  energy  required  to  overcome  headloss  (kW-hr  per  m  of  water;  see  appendix  B) 
Priceeiec  =  price  of  electricity  ($/kW-hr) 

We  can  now  explicitly  write  objective  f>  as  a  function  of  decision  variables: 
Operating  cost  =  Material  cost  +  Pump  cost 

HQ,  cin,p)  =  2*  Q*  t*  Cin  *  p/32  *  1000  L/m3  *  Priced0nor  +  2*  Q*  t*  E*  Priceeiec 
HQ,  Ci„,p)  =  2  *  Q*  t*  (Cin  *  p/32  *  1000  L/m 3  *  Priced0nor  +  E  *  Priceelec ) 
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Finally,  we  must  recognize  the  constraints  that  our  problem  domain  imposes.  The  site 
and  technology  models  impose  lower  and  upper  bounds  on  each  of  the  decision  variables 
which  We  can  designate  as  Qmin ?  Qmaxi  dmin,  dmaxi  (- ' m.nttn-  t  />< it v-  P mi u-  and  P max- 

With  our  objectives,  decision  variables,  and  constraints  identified,  we  can  now 
formulate  a  multi-objective  problem:  search  for  all  vectors  x  =  \Q,  d,  Cin,p]  that  optimize 
the  objective  function  vector  f(x)  =  [fi(x),  f2(x)]  subject  to  the  following  constraints: 

Qmin  <Q<  Qmax'i  real-valued 
d„,in  <d<  dmax;  integer- valued 
Cm, min  —  Cin  —  Cin,  maxi  real-valued 
Pmin  —  P  < Pmax ;  real-valued. 

For  this  problem,  HFTW  placement  is  restricted  to  cells  in  column  13  of  the  finite 
difference  grid  shown  in  Figure  3.1;  furthermore,  the  two  HFTWs  are  required  to  remain 
equidistant  from  the  grid  centerline.  Hence,  when  looking  at  Figure  3.1,  well-spacing  (d) 
is  the  vertical  grid  distance  between  HFTWs.  Note  that  d  must  be  integer  because  it 
measures  the  summation  of  vertical  cells  between  the  HFTWs,  and  individual  cell 
dimensions  are  integer  (3,  5,  or  10  meters).  The  remaining  decision  variables  ( Q ,  Cin,p), 
however,  can  accommodate  better  resolution  than  integer  values,  and  therefore,  are 
designated  as  real. 

The  objective  of  this  problem  is  to  find  sets  of  engineering  parameters  [ Q ,  d,  Cin, 
p]  that  yield  an  optimal  trade-off  relationship  between  technology  performance  (with 
respect  to  perchlorate  removal)  and  cost.  Ideally,  we  would  further  constrain  the  solution 
set  to  values  [ Q ,  d,  Cin,  p]  that  yield  downgradient  ClOf  concentrations  that  fall  below 
some  maximum  level;  the  rationale  for  this  constraint  is  that  regulations  typically 
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prescribe  a  maximum  contaminant  level  downgradient  from  the  source.  However,  the 
ability  of  the  technology  model  to  yield  sufficiently  low  CIO4'  levels  is  largely  uncertain; 
this  additional  constraint  would  risk  over-constraining  the  problem  and  yielding  no 
solutions.  Therefore,  instead  of  formulating  downgradient  CIO4'  concentration  as  an 
objective  or  constraint,  we  simply  monitor  downgradient  CIO4"  concentration  to  gain  an 
understanding  of  how  different  technology  implementations  affect  the  relative  magnitude 
of  this  important  parameter. 

There  are  some  important  considerations  in  selecting  a  suitable  search  algorithm  to 
solve  this  problem: 

•  Objectives  lj  and  f2  have  common  decision  variables  \Q,  Cin,p ],  but  the 
mathematical  relationship  between  f  and  f2  is  unknown.  In  other  words  there  is 
no  analytical  expression  for  the  theoretical  Pareto  Front  (PF*). 

•  An  explicit  relationship  between  fi  and  [ Q ,  d,  Cin,  p]  is  unavailable. 

•  Because  a  discrete  plot  of  the  Pareto  front  is  the  only  representation  achievable, 
the  plot  must  adequately  span  the  extremities  of  the  Pareto  front. 

•  Due  to  the  relatively  long  computation  time  in  evaluating  a  decision  variable  set, 
an  efficient  search  algorithm  is  essential. 

Given  these  considerations,  a  genetic  algorithm  seems  appropriate  for  estimating  PF*. 
As  stated  in  Section  2.5,  GAs  adapt  easily  to  multi-objective  problems  because  they 
operate  on  a  population  of  candidate  solutions,  enabling  the  discovery  of  multiple 
nondominated  points  in  a  single  generation. 


54 


3.5  MULTI-OBJECTIVE  GENETIC  ALGORITHM 


3.5.1  MOGA  specifications 

In  developing  an  MOGA  to  solve  the  proposed  problem,  the  following  MOEA 
characteristics  are  desireable: 

•  Pareto-based  to  generate  discrete,  approximated  Pareto-optimal  solutions. 

•  Real-valued  crossover  and  mutation  operators  to  accommodate  possible  ranges  of 
values  for  the  decision  variables. 

•  Fitness  sharing  to  generate  points  that  adequately  span  the  Pareto  front. 

•  Parallelized  to  expedite  computation  by  allowing  multiple  processors  to 
individually  calculate  fitness  for  a  particular  chromosome.  Parallelization  is  a 
practical  necessity  due  to  computation  intensity  of  computing  the  FORTRAN 
flow-and-transport  model. 

Section  2.6  described  several  multi-objective  GAs.  For  this  thesis  effort,  the  author 
implemented  a  version  of  NPGA  largely  due  to  its  success  in  previously  solving  a  multi¬ 
objective  groundwater  remediation  problem  that  is  somewhat  similar  to  the  problem 
being  considered  here  (Erickson  et  al.,  2001).  Also,  the  author  considered  NPGA  easily 
understandable  and  easiest  to  develop  into  low-level  computer  code.  From  this  point, 
forward,  the  algorithm  will  be  referred  to  as  HK-MOGA,  in  which  the  letters  H  and  K 
refer  to  the  software  developers  (Hendricks  and  Knarr;  see  Acknowledgments). 
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3.5.2  Chromosome  encoding 


A  chromosome  consists  of  genes  representing  the  decision  variables  listed  in 
Section  3.4.  In  essence,  a  chromosome  contains  engineering  parameters  for  a  2-well 
HFTW  system.  Chromosomes  also  have  auxiliary  genes  devoted  to  objective  function 
values,  Pareto-ranking,  and  downgradient  CIO4’  concentration,  but  these  auxiliary  genes 
are  not  involved  in  crossover  or  mutation.  Decision  variables  are  real-valued  except  for 
well-spacing,  which  is  integer-valued. 

3.5.3  MOGA  parameters 

HK-MOGA  relies  on  the  following  user-specified  parameters: 

•  Initial  population  size  (IPopJ):  self-explanatory 

•  Mating  pool  size  (\MP\):  sets  the  quantity  of  chromosomes  that  are  selected  for 
crossover  and  mutation. 

•  Number  of  generations  (N):  determines  how  many  cycles  of  selection,  crossover, 
and  mutation  occur. 

•  Niche  radius  (cjsWe):  dimensionless  number  used  in  fitness  sharing;  establishes  a 
circular  niche  around  each  point  in  non-dimens ionalized  objective  space 

•  Mutation  probability  (pm):  self-explanatory 

3.5.3  Initialization 

HK-MOGA  starts  by  randomly  generating  an  initial  chromosome  population  Popo 
of  size  |  Popo|  in  which  all  chromosomes  comply  with  decision  variable  constraints  listed 
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in  section  3.4.  HK-MOGA  stores  the  created  chromosomes  in  a  set  called  PopCUm  ,  so- 
called  because  it  is  a  cumulative  collection  of  chromosomes;  HK-MOGA  continually 
appends  newly  generated  chromosomes  to  Popcum. 

After  Popcum  is  created,  HK-MOGA  evaluates  objective  function  values  for  each 
chromosome.  HK-MOGA  next  uses  the  objective  function  values  to  Pareto-rank  each 
chromosome  in  PopCUm-  A  particular  chromosome  x  receives  a  rank  that  equals  the 
number  of  chromosomes  that  dominate  chromosome  x.  Hence,  all  nondominated 
solutions  have  rank  =  0,  and  low  rank  corresponds  to  high  fitness.  This  Pareto-ranking 
method  mimics  the  method  of  Fonseca  and  Fleming  (1993). 

3.5.4  Selection 

Once  chromosomes  are  Pareto-ranked,  HK-MOGA  starts  the  generation 
subroutine,  which  begins  by  selecting  chromosomes  to  be  copied  to  a  reservoir  called  a 
mating  pool  ( MP )  where  they  await  crossover  and  mutation.  The  user-specified 
parameter  \MP\  limits  the  number  of  chromosomes  in  MP. 

Two  properties  that  drive  selection  in  HK-MOGA  are  1)  Pareto  rank  and  2) 
“crowding”  in  objective  space.  As  long  as  space  in  MP  is  sufficient,  HK-MOGA  copies 
all  rank-zero  chromosomes  to  MP,  then  all  rank-l’s,  then  all  rank-2 ’s,  etc.  Simply  put, 
when  space  is  sufficient,  Pareto  rank  drives  selection.  This  method  ensures  that  better, 
low-ranking  chromosomes  receive  priority  for  selection. 

However,  as  low-ranking  chromosomes  are  progressively  copied,  space  in  MP 
depletes.  Eventually  HK-MOGA  encounters  chromosomes  of  some  rank  k  whose 
quantity  exceeds  remaining  space.  Because  all  /c-rank  chromosomes  have  equal  rank,  and 
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therefore  fitness,  HK-MOGA  resolves  the  dilemma  of  which  chromosomes  to  select  by 
picking  chromosomes  that  are  least  crowded  in  objective  space.  The  purpose  of  this 
selection  strategy  is  to  develop  chromosomes  that  map  to  relatively  uninhabited  sections 
of  the  Pareto  front,  which  for  the  decision  maker  means  more  diverse  tradeoff  options. 
The  metric  for  assessing  crowding  in  objective  space  is  the  niche  count.  This  technique 
is  called  equivalence  class  sharing,  which  was  originally  described  by  Horn  et  al. 
(1994). 

Let  x;  be  some  A:- rank  chromosome  such  that  i  =  1,2,3...  [number  of  rank  k 
chromosomes],  and  assume  that  the  number  of  A: -rank  chromosomes  exceeds  remaining 
space  in  the  mating  pool.  Also,  let  chromosome  x,  be  any  chromosome  in  Popcum  where  j 
=  1,2,3...  [size  of  Popcum] .  HK-MOGA  searches  Popcum  for  the  most  current  maximum 
and  minimum  values  of  both  objective  functions  fi  and  f2.  These  maximum  and 
minimum  values  are  subsequently  used  to  nonnalize  objective  function  values  for  every 
X;  e  {rank  k  chromosomes}  and  Xj  e  Popcum  as  follows: 

ffi  =  (fii  -  fimin)  /  (fimax  “  fimin) 
ffi  =  (fli  -  fimin)  /  (f2max  “  fimin) 
ffj  =  (flj  “  fimin)  /  (fimax  “  fimin) 
ffj  =  (f2j  “  fimin)  /  (fimax  “  fimin) 

where 

fi'i  =  dimensionless  value  of  fi  based  on  x;e  {rank  k  chromosomes} 
fi'i  =  dimensionless  value  of  f>  based  on  x;e  {rank  k  chromosomes} 
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fi'j  =  dimensionless  value  of  fi  based  on  XjePopcum 

fi'j  =  dimensionless  value  of  f>  based  on  XjePopcum 

fa  =  value  of  fj  based  on  xje  {rank  k  chromosomes} 

f2i  =  value  of  fj  based  on  x;e  {rank  k  chromosomes} 

fij  =  value  of  fi  based  on  XjePopcum 

f2j  =  value  of  fz  based  on  XjePopcum 

fi  min  =  minimum  value  of  fi  within  Popcum 

fimax  =  maximum  value  of  fi  within  Popcum 

f2min  =  minimum  value  of  f2  within  PopCUm 

f2max  =  maximum  value  of  f>  within  Popcum 

This  normalization  makes  both  objective  function  values  dimensionless,  which  is  helpful 
due  to  the  incommensurable  units  of  both  objective  functions  (mass  and  dollars).  The 
distance  d,,  between  points  (ffo,  fz'i)  and  (fi'j,  f^'j)  in  dimensionless  objective  space  is 
calculated  as 

dij  =  [  (fi'r-  fi))2  +  (f2'i  -  f2))2  ],/2 

Distance  d;j  and  the  niche  radius  ashare  are  then  used  to  compute  the  sharing 
function  (equation  2.11): 

Sh(dij)  —  1  —  dij  /cishare  for  dij  £  CTshare 

=  0  for  dij  >  CJ share 

The  parameter  ashare  basically  defines  the  radius  of  a  circle  around  point  (ffo,  IV,);  points 
inside  the  circle  contribute  to  crowding,  and  points  outside  the  circle  do  not.  Sh(dq)  is  a 
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metric  for  assessing  the  proximity  or  crowding  between  point  (ffj,  ffj)  and  some  other 
point  (ffj,  ffj).  If  (ffj,  ffj)  lies  within  the  circle  surrounding  (f',,  (V,)  (i.e.  d,,  <  ashare), 
then  the  sharing  function  assumes  a  value  such  that  0  <  Sh(dy)  <  1 .  The  closer  the  two 
points  are  to  each  other,  the  higher  the  value  of  the  sharing  function;  the  maximum  value 
Sh(dij)  =  1  indicates  the  two  points  overlap.  If  point  (ffj,  ffj)  lies  outside  the  circle 
surrounding  (ff ;,  f2'i)  (i.e.  d,,  >  ashare),  then  crowding  is  negligible  (Sh(dij)  =  0). 

The  niche  count  m,  for  chromosome  x,e  {rank  k  chromosomes}  is  computed 
according  to  equation  2.12: 

/U;  ^jxjePopcum  Sh(dij) 

Because  i?n  is  a  summation  of  sharing  function  values,  it  provides  an  overall  measure  of 
how  “crowded”  chromosome  x;  is  in  objective  space.  A  high  niche  count  implies  a  high 
degree  of  crowding,  and  vice  versa  (see  Figure  4).  After  calculating  the  niche  count  m, 
for  each  x,e  {rank  k  chromosomes},  HK-MOGA  fills  remaining  slots  in  MP  with 
chromosomes  having  the  lowest  m,  values.  HK-MOGA  preferentially  selects 
chromosomes  with  low  niche  counts  to  improve  chances  of  generating  points  in  less- 
occupied  regions  of  the  Pareto  front. 
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Figure  4.  Example  of  niching  strategy.  Point  B  is  less  crowded  than  point  A. 

3.5.5  Crossover  and  Mutation 

Chromosomes  in  MP  proceed  to  crossover  and  mutation.  HK-MOGA  relies  on 
several  crossover/mutation  operators,  which  Garrett  (1999)  coded  for  his  thesis 
investigation.  Garrett’s  (1999)  computer  code  was  incorporated  into  HK-MOGA  code. 
Although  Section  2.6.4  addresses  mating  restriction  during  crossover,  mating  restriction 
does  not  appear  to  be  a  critical  component  of  MOGAs  investigated  in  the  literature 
review,  and,  as  stated  before,  no  sound  theory  justifies  its  inclusion  in  MOGAs 
(Veldhuizen  and  Lamont,  2000).  Therefore,  the  author  did  not  incorporate  mating 
restriction  with  the  crossover. 

All  chromosomes  in  MP  participate  in  crossover.  That  is,  for  i  =  1  to  \MP\, 
chromosome  x,  crosses  over  with  xr,  where  xr  is  randomly  chosen  from  MP  (pool 
crossover  operates  differently).  HK-MOGA  employs  the  following  crossover  operators: 
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Whole  arithmetical  crossover:  linearly  combines  all  corresponding  genes  of  x;  and 
xr,  to  create  new  chromsomes  x'i  and  x'2;  crossover  applies  to  the 
“whole”  chromosome  (i.e.  all  genes  of  x;  and  xr)  (Michaelwicz, 

1996:112,  128).  HK-MOGA  randomly  retains  x'i  and  discards  x'2. 

Simple  crossover:  randomly  selects  a  gene  and  swaps  it  between  Xj  and  xr  to  make 
x'i  and  x'2  (Michaelwicz,  1996: 112).  HK-MOGA  randomly  retains 
xhand  discards  x'2. 

Heuristic  crossover:  uses  chromosomes  x;  and  xr  to  make  a  single  offspring  x'i 

such  that  x'i  =  R(  xr  -  \,)  +  xr.  The  value  R  is  a  unifonn  random  number 
between  0  and  1 ,  and  the  rank  of  xr  is  the  same  or  less  than  the  rank  of  x, 
(Michaelwicz,  1996:112). 

Pool  crossover:  randomly  copies  alleles  from  chromosomes  in  MP  and  assembles 
the  alleles  to  make  x'i. 

HK-MOGA  selects  a  particular  crossover  operators  based  upon  an  adaptive  probability 
distribution  (Garrett,  1999).  At  the  first  generation,  all  crossover  operators  have  equal 
probability  of  selection.  For  all  following  generations,  the  selection  probability  for  a 
particular  operator  “adapts”  or  adjusts  based  upon  the  attributes  of  the  new  chromosome 
x'i.  If  x'i  dominates  x,,  then  the  crossover  operator  was  successful  in  increasing  fitness, 
and  its  selection  probability  consequently  increases  in  the  next  generation.  Conversely,  if 
x;  dominates  x'i,  then  crossover  was  unsuccessful,  and  its  probability  of  selection 
decreases.  If  neither  chromosome  dominates  the  other,  the  operator’s  selection 
probability  stays  the  same. 
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Crossover  creates  new  chromosomes  that  are  then  susceptible  to  mutation. 
Mutation  is  controlled  by  the  user-specified  mutation  probability  (pm).  For  each  new 
chromosome,  HK-MOGA  selects  a  random  number  (r.  0  <  r  <  1)  from  a  uniform 
distribution.  If  r  <  pm,  then  one  of  3  mutation  operators  affects  the  new  chromosome; 
otherwise,  mutation  does  not  occur.  HK-MOGA  randomly  selects  which  mutation 
operator  based  on  the  same  adaptive  probability  distribution  described  previously: 

Uniform  mutation:  resets  a  particular  gene  to  a  random  value  between  specified 
maximum  and  minimum  values  (Michaelwicz,  1996: 111,  127). 

Boundary  mutation:  resets  a  particular  gene  to  either  its  specified  maximum  or 
minimum  value  (Michaelwicz,  1996:127-128). 

Non-uniform  mutation:  modifies  a  particular  gene  by  some  random  value  whose 
magnitude  decreases  probabilistically  towards  0  as  the  current 
generation  number  approaches  the  maximum  generation  number 
(Michaelwicz,  1996:111,  128). 

Crossover  and  mutation  ultimately  create  a  new  chromosome  population  Popnew,  whose 
size  equals  the  mating  pool  size  \MP\. 

3.5.6  Evaluation  and  Pareto  Ranking 

HK-MOGA  evaluates  all  members  of  Popnew  and  appends  them  to  PopCUm-  Thus, 
Pop  cum  keeps  all  chromosomes  from  past  generations  and  inherits  new  ones.  We  see  then 
that  Popcum  is  a  set  of  accumulated  chromosomes,  and  its  cardinality  is  prescribed  by  the 
following  formula: 

|Pop  cum  |  =  N-\MP\  +  |Pop0| 
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HK-MOGA  ranks  each  chromosome  in  PopCUm  as  described  previously,  and  the 
generation  cycle  restarts  with  selection. 


3.6  OPTIMIZATION  PROCEDURE 


The  author  designed  optimization  runs  to  achieve  the  second  and  third  objectives 
mentioned  in  Section  1.2.  Four  separate  optimization  runs  were  performed  using  HK- 
MOGA  described  in  section  3.5  to  produce  an  estimated  Pareto  set  (PknOWn)  and  Pareto 
front  (PFknown)  for  different  time  spans  and  site  conditions.  Runs  1  and  3  simulate 
treatment  periods  of  300  days  and  600  days,  respectively,  using  site  data  for  Site  4, 
Nevada  (Parr,  2002): 

Aquifer  Characteristics 

Hydraulic  Conductivity  =  7.6  m/day  Hydraulic  Gradient  =  0.01 

Source  Characteristics 

Oxygen  Concentration  =  2.8  mg/L  Nitrate  Concentration  =  60  mg/L 

Perchlorate  Concentration  =  330  mg/L 


Runs  2  and  4  simulate  treatment  periods  of  300  days  and  600  days,  respectively,  using 
site  data  for  Site  2,  California  (Parr,  2002): 

Aquifer  Characteristics 

Hydraulic  Conductivity  =  2.59  m/day  Hydraulic  Gradient  =  0.001 

Source  Characteristics 

Oxygen  Concentration  =  0.55  mg/L  Nitrate  Concentration  =  0.5  mg/L 

Perchlorate  Concentration  =160  mg/L 


As  mentioned  in  Section  1.2,  one  of  the  research  objectives  was  to  perform  the 
optimization  under  “various  contaminated-site  conditions.”  Key  parameters  that  establish 
different  conditions  between  the  two  sites  are  hydraulic  conductivity,  regional  hydraulic 
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gradient  and  initial  source  concentration.  Site  4  has  approximately  triple  the  hydraulic 
conductivity,  ten  times  the  hydraulic  gradient,  and  double  the  source  concentration  of  Site 
2;  Site  4  also  has  larger  initial  concentrations  of  competing  electron  acceptors  (O2  and 
NO3")  than  Site  2. 

The  following  parameters  are  used  in  all  four  runs: 

Decision  variable  constraints 

Qmin  =  10  m3/day  Qmax  =  150  m3/day  dmin  =  3  m  dmax  =  165  m 

Cin,min  0  Hlg/L  Cin,  max  1,000  Hlg/L 

Pmin  0  Pmax  =  32 

Cost  coefficients 

Pricedonor  =  $2.666x  10"6  per  mg  donor  Priceeiec  =  $0,067  per  kW-hr 

HK-MOGA  parameters 

|Pop0|  =  50  \MP\  =  10  N=  100  ashare  =  0.4  pm  =  0.01 

Aquifer  parameters 

Porosity  =  0.30  Retardation  factor  for  acetate  (CH3COO')  =  1 .48 

Kinetic  parameters  (see  Appendix  A) 

kmax  =  0.2 1  mg  donor/mg  biomass/day  Ksdon  =  10.0  mg/L 

Ksoxy  =10.0  mg/L  Ksmt  =15.0  mg/L  Ksper=  20.0  mg/L 

K°xy  =10.0  mg/L  K,nit  =15.0  mg/L 

Ybiomass  =  0.25  mg  biomass/mg  donor 

Foxy  =  0.83  mg  oxygen/mg  donor  Fnit  =1.3  mg  nitrate/mg  donor 

Fper  =  1.45  mg  oxygen/mg  donor  b  =  0.002  day'1  Xm;n  =  0.01  mg/L 

The  Qmax  was  selected  as  an  appropriate  real-world  value  that  would  not  result  in 

excessive  drawdown  at  the  sites.  Values  for  dmin  and  dmax  are  based  upon  the  dimensions 

of  the  site  model.  The  value  for  Pricedonor  is  based  upon  an  estimated  bulk  cost  of 

$286.20  per  55-gallons  of  a  50/50  mixture  of  acetic  acid  (CH3COOH)  and  water.  The 

value  for  Priceeiec  came  from  a  U.S  Department  of  Energy  (USDOE)  website  that  lists 

average  electricity  prices  for  commercial  consumers  in  Nevada  in  the  year  2000 

(USDOE,  2003).  All  values  for  the  kinetic  parameters  are  the  same  as  the  baseline  values 

that  Parr  used  (2002)  except  for  the  biomass  decay  rate,  which  was  changed  to  a  value  of 

b  =  0.002  day'1  (originally  b  =  0.01  day'1)  to  simulate  less  biomass  decay.  HK-MOGA 
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parameters  values  are  based  on  both  knowledge  of  similar  parameter  values  seen  in  the 
literature  and  off-line  experiments  with  the  HK-MOGA  software. 

The  HK-MOGA  software  performs  the  evaluation  function  (Section  3.5.6)  by 
executing  the  FORTRAN-coded  technology  model  and  computing  the  cost  formula  for 
each  new  chromosome.  Execution  of  the  technology  model  is  the  most  computationally 
expensive  and  time-intensive  activity  in  the  program.  HK-MOGA,  which  is  coded  in 
C++,  uses  Message-Passing  Interface  (MPI)  to  enable  parallel  computation  among 
several  Aspen  dual-processor  machines.  Each  Aspen  machine  has  1-GB  memory  and 
two  1-GHz  Pentium  III  processors  that  can  separately  evaluate  the  technology  model  for 
each  chromosome.  Evaluations  for  300-day  and  600-day  scenarios  typically  lasted  about 
8  and  17  minutes,  respectively,  on  each  processor.  The  program  was  executed  with 
Redhat  LINUX  version  7.3  and  MPI  version  1.2. 7.1. 

In  addition  to  decision  variable  and  objective  function  values,  the  HK-MOGA 
software  also  outputs  the  maximum  CIO4'  concentration  among  all  cell  layers  of  column 
20  of  the  finite  difference  grid  (the  monitoring  wells  shown  in  Figure  3).  As  mentioned 
before,  the  purpose  of  including  this  output  is  to  assess  the  ability  of  the  treatment 
technology  to  meet  regulatory  requirements. 
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4.0  RESULTS  AND  ANALYSIS 


4.1  INTRODUCTION 

In  this  chapter  we  present  and  discuss  the  results  obtained  from  using  HK-MOGA 
described  in  Section  3.5  to  optimize  both  perchlorate-mass  destruction  and  operating  cost 
incurred  based  upon  simulations  of  Parr’s  (2002)  in  situ  bioremediation  technology 
model.  First  we  compare  Runs  1  and  2  (300-day  scenarios)  for  Sites  4  (in  NV)  and  2  (in 
CA)  with  respect  to  their  estimated  Pareto  fronts  (PFkn0Wn)  and  their  downgradient 
perchlorate  concentration  ([CIO4"])  data.  Next,  a  similar  analysis  is  presented  for  Runs  3 
(Site  4)  and  4  (Site  2)  for  600  days  of  technology  operation. 


4.2  RESULTS  FOR  300-DAY  TECHNOLOGY  OPERATION 

Figures  5  and  6  show  PF known  for  Runs  1  and  2,  respectively,  at  generation  N  = 
100.  Each  graph  displays  evaluations  of  CIO4’  mass  removal  (fi)  and  operating  cost  (13) 
assuming  a  300-day  time  span.  At  N  =  100,  Runs  1  and  2  generated  429  and  461 
nondominated  points,  respectively,  out  of  a  total  of  1,050  points.  For  both  runs,  HK- 
MOGA  discovered  the  same  upper  and  lower  extremities  with  respect  to  operating  cost, 
and  for  both  runs  PFknown  acquired  sufficient  definition  (i.e.  converged)  by  N  ~  50.  The 
niching  strategy  described  in  Section  3.5.4,  with  ashare  =  0.4,  was  effective  in  generating  a 
representative  span  of  points  between  the  nondominated  boundary  points. 


67 


f2:  Operating  cost 
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Figure  5.  Run  1  est.  Pareto  front  (t  =300  days;  Site  4  parameters). 


fl:  Perchlorate  mass  removed  (kg) 

Figure  6.  Run  2  est.  Pareto  front  (t  =300  days;  Site  2  parameters). 
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Both  nondominated  fronts  for  Runs  1  and  2  exhibit  a  nonlinear  relationship  between  fi 
and  f>,  with  dfVdfi  increasing  with  increasing  values  of  f.  This  indicates,  as  might  be 
anticipated,  that  increments  of  perchlorate  mass  removal  become  increasingly  expensive. 
Section  3.4  indicated  that  operating  cost  consisted  of  two  parts:  the  cost  of  electron  donor 
and  the  cost  of  operating  the  pumps.  If  Q  assumes  the  value  Qmax,  then  the  cost  of 
operating  the  pumps  is  maximized  at  about  $780  for  the  300-day  period,  which  is  a 
negligible  component  of  total  costs,  which  can  be  seen  from  Figures  5  and  6  to  be  in  the 
tens  of  thousands  of  dollars  for  any  mass  removal  much  greater  than  zero.  Therefore, 
high  operating  costs  are  essentially  attributed  to  usage  of  electron  donor  alone. 

One  explanation  for  the  increasing  incremental  costs  of  mass  removal  is  that  at 
high  mass  removal  rates  (which  are  needed  to  achieve  high  total  mass  removal)  the 
perchlorate  reduction  reaction  becomes  limited,  either  by  kinetics  or  by  biomass.  Thus, 
increasing  donor  is  needed  to  achieve  equivalent  rates  of  perchlorate  reduction,  as 
perchlorate  concentration  increases.  In  the  limit,  the  rate  of  perchlorate  reduction  is 
maximized  and  additional  donor  addition  (resulting  in  increased  cost)  has  no  effect  on  the 
rate  or  extent  of  perchlorate  reduction. 

Runs  1  and  2  were  for  two  different  sites,  with  different  hydraulic  conductivities 
and  initial  CKV  concentrations.  Figure  7  plots  PFknown  for  both  Runs  1  and  2  on  the 
same  set  of  axes.  Although  PF|<nown  for  Runs  1  and  2  share  the  same  range  of  operating 
cost  values,  nondominated  points  for  Run  2  indicate  considerably  less  mass-removal  than 
those  for  Run  1.  In  fact,  Figure  7  shows  that,  for  the  same  operating  cost  of  $240,771, 
Run  1  achieved  a  maximum  mass-removal  value  of  5,170  kg  of  CKV  destroyed 
compared  to  only  175  kg  for  Run  2.  The  large  disparity  in  CKV  removal  performance  is 
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almost  certainly  attributable  to  differences  in  hydrogeological  properties.  Hydraulic 
conductivity,  which  quantifies  how  easily  water  flows  through  the  soil  for  a  given 
hydraulic  gradient,  is  about  3  times  smaller  in  Site  2  than  Site  4  (2.59  vs.  7.60  m/day); 
hydraulic  gradient  is  also  10  times  lower  for  Site  2  (0.001  vs.  0.010).  The  combination  of 
these  parameters  results  in  a  regional  groundwater  Darcy  velocity  for  Site  2  that  is  ~  1/30 
that  of  Site  4.  This  low  Darcy  velocity,  in  turn,  causes  higher  interflow/recirculation 
between  the  two  HFTWs  for  Site  2.  Increased  recirculation  means  that  for  a  given 
pumping  rate,  there  is  less  capture  of  upgradient  contaminant  by  the  treatment  system, 
and  less  mass  removal.  Another  potential  reason  for  the  mass-removal  disparity  between 
the  two  sites  is  the  difference  in  initial  CIO4'  concentration:  330.0  mg/L  for  Run  1  (Site  4) 
versus  160.0  mg/L  for  Run  2  (Site  2).  Clearly,  contaminant  mass  removal  is  likely  to  be 
higher  if  higher  concentrations  of  contaminant  are  present  at  the  site. 
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The  technology  model  also  provided  information  about  downgradient 
concentration  ([CIO4'])  by  comparing  concentration  values  in  all  cells  of  column  20  of 
the  discretized  site  model  (Figure  3)  and  then  outputting  the  maximum  value.  The 
purpose  of  observing  this  attribute  was  to  assess  the  bioremediation  technology’s 
potential  to  attain  hypothetical  regulatory  limits  on  downgradient  contaminant  levels. 

Figure  8  shows  the  plot  of  maximum  downgradient  concentration  versus  CIO4' 
mass-removal  for  Run  1 .  The  following  observations  can  be  made  regarding  this  figure: 

•  In  general,  solutions  (both  nondominated  and  dominated)  with  higher  mass 
removal  tend  to  have  lower  concentration  measurements,  as  would  be  expected. 
Note,  however,  that  there  are  large  individual  variations  within  this  general  trend. 

•  Nondominated  points  with  the  lowest  downgradient  [CKV]  had  very  similar 
genotypes  ( i.e .  decision  variable  values);  Q,  Cin,  and p  assumed  values  at  or  near 
Qmax,  Cin, max,  and  pmax,  respectively,  indicating  large  mass  per  time  of  electron 
donor  injected.  Consequently,  the  objective  function  values  associated  with  these 
points  -tended  to  be  relatively  high:  mass  removal  «  5,170  kg  and  operating  cost « 
$241,000.  Well  spacing  (d)  was  typically  25  meters.  Downgradient  [CKV]  for 
these  particular  nondominated  solutions  varied  from  26  to  30  mg/L;  greatly 
exceeding  the  state-specified  action  levels  in  the  pg/L-range  that  were  discussed 
in  Section  2.1. 

•  Many  dominated  solutions,  despite  being  inferior  in  terms  of  mass  removal  and 
cost,  yielded  lower  CIO4'  downgradient  levels  than  their  nondominated 
counterparts.  Interestingly,  the  lowest  concentration  for  Run  1  (24.1  mg/L) 
originates  from  a  solution  that  is  inferior  to  326  other  solutions. 


71 


Figure  9  shows  the  plot  of  maximum  downgradient  concentration  versus  CIO4'  mass- 
removal  for  Run  2.  The  following  observations  apply  to  Figure  9: 

•  The  plot  of  all  points  (both  nondominated  and  dominated)  reveals  no  discernible 
relationship  between  [CIO4']  and  mass  removal. 

•  A  chain  of  nondominated  points  with  mass  removal  >  40.0  kg  had  relatively 
narrow  range  of  downgradient  [CIO4’]  values  (8.4  -  8.7  mg/L);  genotypes 
associated  with  [CKV]  «  8.4  mg/L  were  nearly  homogeneous;  Q,  Cin,  and  p  had 
values  at  or  near  Qmax,  Cinymax,  and  pmax,  respectively,  again  indicating  that  mass  of 
electron  donor  injected  per  time  was  large.  Well  separation  (d)  was  typically  25 
meters. 

•  Ironically,  certain  nondominated  points  associated  with  nearly-zero  mass  removal 
had  lower  [CIO4']  values  than  other  nondominated  points  with  higher  removals  of 
CIO4'  mass.  These  points  had  minimum  values  for  Q  as  well  as  Cm  and/or  p 
minimized,  indicating  relatively  low  rates  of  donor  mass  injected  per  time.  Well 
spacing  ( d)  varied  from  35  -  45  meters,  and  downgradient  [CKV]  was  «  3xl0'4 
mg/L,  which  complies  with  state-specified  action  levels  mentioned  Section  2.1. 
These  low  concentration  measurements  may  be  attributed  to  interflow  (Section 
2.2),  which  is  a  function  of  the  design  parameters  ( Q  and  d)  and  the  regional 
groundwater  Darcy  velocity.  Site  2  hydraulic  conductivity  is  ~  1/3  that  of  Site  4, 
and  Site  2  hydraulic  gradient  is  1/10  that  of  Site  4;  this  combination  results  in  a 
regional  groundwater  Darcy  velocity  for  Site  2  that  is  30  times  smaller  than  that 
of  Site  4.  Although  the  low  Q  for  these  solutions  would  tend  to  reduce  interflow, 
the  combination  of  low  Darcy  velocity  and  small  well  spacing  ultimately 
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Max.  CI04' concentration  (mg/L) 


increases  interflow.  High  interflow  means  that  although  very  little  contaminant 
mass  is  treated,  the  contaminated  water  that  is  treated  ends  up  passing  through  the 
treatment  wells  multiple  times,  and  downgradient  contaminant  concentrations 
would  be  low. 

•  Figure  9  shows  a  large  number  of  dominated  points  that  yielded  lower  C104‘ 
downgradient  concentrations  than  other  nondominated  points.  In  fact, 
nondominated  points  with  mass  removal  >  40.0  kg  seem  to  establish  an  upper 
bound  on  measured  concentration,  with  many  dominated  points  lying  below. 


f  1 :  Perchlorate  mass  removed  (kg) 

Figure  8.  Run  1:  maximum  downgradient  perchlorate  cone.  ([C104~])  vs.  C104'  mass  removed. 
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fl:  Perchlorate  mass  removed  (kg) 

Figure  9.  Run  2:  maximum  downgradient  perchlorate  cone.  ([CIO4])  vs.  CIO4"  mass  removed. 


4.3  RESULTS  FOR  600-DAY  TECHNOLOGY  OPERATION 

Figures  10  and  1 1  show  the  plots  of  PFkn0wn  for  Runs  3  and  4,  respectively,  after 
100  generations.  These  plots  are  very  similar  in  shape  to  the  respective  300-day  runs 
(Figures  5  and  6),  which  is  logical  due  to  identical  optimization  parameters  (except  time). 
Upper  boundary  values  for  mass  removal  and  operating  cost  are  larger,  of  course,  due  to 
the  longer  time  span.  The  same  analysis  used  in  section  4.1  can  describe  the  nonlinear 
increase  in  incremental  costs  with  mass  removal  that  is  seen  here. 

Figure  12  plots  PFknown  for  Runs  3  and  4  on  the  same  set  of  axes.  Figure  12 
reveals  the  same  relationship  that  Figure  7  does:  Site  4  parameters  yield  higher  CIO4' 
mass  removal  than  Site  2  parameters  for  the  same  range  of  operating  cost  for  reasons 
mentioned  previously  with  the  300-day  simulations  (different  hydrogeolocial  parameters 
and  initial  concentration). 
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f  1 :  Perchlorate  mass  removed  (kg) 

Figure  10.  Run  3  est.  Pareto  front  (t  =600  days;  Site  4  parameters). 
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fl:  Perchlorate  mass  removed  (kg) 

Figure  11.  Run  4  est.  Pareto  front  (t  =600  days;  Site  2  parameters). 
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Figure  12.  Runs  3  and  4  est.  Pareto  fronts  (t  =600  days). 


Figure  13  plots  maximum  downgradient  concentration  versus  CIO4'  mass  removal  for 
Run  3  after  100  generations.  The  following  observations  can  be  made  regarding  Figure 
13: 

•  Just  like  Run  1 ,  solutions  (both  nondominated  and  dominated)  with  higher  mass 
removal  generally  provide  lower  concentration  measurements.  Nondominated 
solutions  appear  to  follow  this  trend  better  than  the  dominated  solutions. 

•  660  solutions  had  measured  downgradient  [C104‘]  that  was  less  than  the  lowest 
value  reported  for  Run  1  (24.1  mg/L);  283  of  these  solutions  were  nondominated. 
This  seems  to  indicate  that,  for  Site  4  parameters,  downgradient  CIO4'  levels 
improve  for  longer  treatment  periods. 

•  Nondominated  points  with  lowest  downgradient  [CIO4']  had  very  similar 
genotypes;  Q,  Cin,  and  p  values  at  or  near  Qmax,  Cin,max,  and  pmax,  respectively, 
indicating  large  mass  per  time  of  electron  donor  injected  and  high  operating  cost 
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(«$48 1,500).  For  these  particular  points,  well  spacing  (d)  was  uniformly  25 
meters,  and  downgradient  [CIO4’]  was  «0.154  mg/L.  Despite  being  an 
improvement  over  Run  1 ,  this  level  still  exceeds  the  state-specified  action  levels 
discussed  in  Section  2.1,  which  are  in  the  pg/L  range. 

•  Numerous  dominated  points,  despite  being  inferior  in  terms  of  mass  removal  and 
cost,  outperformed  other  nondominated  points  in  achieving  lower  downgradient 
[CIO4'].  The  lowest  concentration  for  Run  3  (0.022  mg/L)  originates  from  a 
solution  inferior  to  330  other  solutions  in  terms  of  mass  removal  and  cost. 


fl:  Perchlorate  mass  removed  (kg) 

Figure  13.  Run  3:  maximum  downgradient  perchlorate  cone.  (|C104"])  vs.  C104‘  mass  removed. 
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Figure  14.  Run  4:  maximum  downgradient  perchlorate  cone.  ([C104])  vs.  C104  mass  removed. 


Figure  14  plots  maximum  downgradient  concentration  versus  CKV  mass-removal  for 
Run  2.  The  following  observations  apply  to  Figure  14: 

•  The  plot  of  all  points  (both  nondominated  and  dominated)  reveals  no  discernible 
relationship  between  [CKV]  and  mass  removal. 

•  Nondominated  solutions  in  the  mass-removal  range  of  >1,000  kg  provided  lower 
[CKV]  measurements  as  mass  removal  increased. 

•  Just  like  Run  2,  Run  4  shows  that  certain  nondominated  points  with  negligible 
mass  removal  showed  lower  [CKV]  than  other  nondominated  points  with  higher 
mass  removal.  These  points  had  low  Q  values  and  one  or  both  of  C„,  and  p 
minimized,  indicating  small  mass  per  time  of  electron  donor  injected.  Well 

■j 

spacing  (d)  was  minimized  (3  meters),  and  downgradient  [CKV]  was  «  9x1  O' 
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mg/L,  which  falls  in  the  pg/L  range  of  state-specified  action  levels  discussed  in 
Section  2.1.  These  low  concentrations  can  also  be  attributed  to  high  interflow 
(due  to  the  low  hydraulic  conductivity  and  hydraulic  gradient  at  the  site)  as 
previously  explained  for  Run  2. 

•  A  multitude  of  dominated  points  yielded  better  ClOfi  downgradient 

concentrations  than  other  nondominated  points.  In  fact,  nondominated  points  in 
the  mass-removal  range  of  >1,000  kg  seem  to  establish  an  upper  bound  on 
measured  concentration,  with  many  dominated  points  lying  below. 


4.4  RESULTS  SUMMARY 

Runs  1,  2,  3,  and  4  were  performed  to  fulfill  the  objectives  listed  in  Section  1.2.  For  each 
of  the  four  runs,  HK-MOGA  estimated  a  Pareto  front  (Figures  5,  6,  10,  and  1 1)  that 
shows  a  nonlinear,  increasing  relationship  between  mass  removal  and  operating  cost. 
Figures  7  and  12  demonstrate  that  the  bio  remediation  technology  was  less  effective  at 
removing  Cl  Of  mass  under  Site  2  conditions  than  under  Site  4  conditions.  For  Site  4 
conditions,  increasing  the  remediation  time  span  from  300  to  600  days  enabled  the 
treatment  technology  to  attain  lower  downgradient  ClOfi  concentrations;  under  Site  2 
conditions,  however,  increasing  the  time  span  rendered  no  discernible  improvement  with 
respect  to  downgradient  concentration.  High  mass  removal  of  ClOfi  is  not  necessarily 
coincident  with  lower  ClOfi  concentration,  especially  under  Site  2  conditions. 
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5.0  CONCLUSIONS 


5.1  SUMMARY 

In  this  thesis  effort,  we  developed  and  presented  a  method  for  optimizing  applications  of 
an  innovative  in  situ  bioremediation  technology  that  uses  HFTWs  to  remediate 
perchlorate-  contaminated  groundwater.  The  method  involves  coupling  a  multi-objective 
genetic  algorithm  with  the  technology  model  developed  by  Parr  (2002).  To  meet  the 
research  goal  and  achieve  the  research  objectives,  it  was  necessary  to  use  the  technology 
model  with  appropriate  parameters  to  simulate  applications  of  the  technology  at 
representative  sites,  formulate  a  multi-objective  problem,  develop  a  multi-objective 
optimization  algorithm  that  would  be  suitable  to  solve  the  problem  by  coupling  with  the 
technology  model,  and  select  appropriate  HK-MOGA  parameters.  The  method  was 
applied  to  two  sites  having  different  hydrogeological  conditions.  HK-MOGA  yielded 
system  designs  (flow  rate,  well  spacing,  injected  nutrient  concentration,  and  injection 
pulse  duration)  that  provide  valuable  insights  into  the  tradeoffs  between  1)  the 
technology’s  performance,  defined  by  contaminant  mass  removal  and  2)  the  operating 
costs  incurred  from  implementing  the  technology. 

The  software  used  in  this  study  was  an  adaptation  of  Garrett’s  (1999)  GA  code.  In 
combination  with  Parr’s  (2002)  technology  model,  HK-MOGA  software  can  serve  as  a 
tool  for  both  optimizing  future  applications  of  this  innovative  bioremediation  technology, 
and  helping  us  to  better  understand  how  HFTWs  can  be  used  in  conjunction  with  in  situ 
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biodegradation.  This  study  contributes  to  efforts  taken  to  resolve  groundwater 
contamination  problems  caused  by  CIO4"  releases  across  the  United  States. 


5.2  CONCLUSIONS 

•  The  multi-objective  genetic  algorithm  (MOGA)  developed  in  this  study  appeared 
useful  for  detennining  technology  design  parameters  for  in  situ  bioremediation  of 
perchlorate  using  HFTWs  to  minimize  cost  (defined  as  operating  cost)  and 
maximize  technology  performance  (defined  as  perchlorate  mass  removal).  HK- 
MOGA  determined  various  sets  of  design  parameters  ( Q ,  d,  Cin,  p)  that  provided 
a  decision  maker  with  combinations  of  cost  and  mass  removal  that  were  “Pareto 
optimal”,  that  is  to  say,  nondominated  by  other  potential  solutions.  This  set  of 
solutions  allows  a  decision  maker  to  select  a  system  design  based  on  the 
weighting  of  the  relative  importance  of  performance  versus  cost.  One 
disadvantage  of  HK-MOGA  is  the  need  for  relatively  extensive  computer 
resources  (time  and  CPU  power)  to  evaluate  the  technology  model.  Also,  the 
selection  of  HK-MOGA  parameters  (|Popo|,  \MP\,  N,  as|iarc,  pm)  was  arbitrary,  and 
based  only  on  the  experience  and  judgment  of  the  investigator.  It  is  unclear 
whether  selection  of  a  different  set  of  parameters  would  result  in  a  different 
outcome. 

•  Pareto  fronts  generated  for  all  four  simulations  indicated  that  the  incremental 
operating  cost  of  the  technology  increases  as  the  technology  removes  more 
perchlorate  mass. 
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•  For  a  contaminated  site  with  relatively  high  conductivity  and  regional  gradient, 
the  choice  of  technology  operating  time  affected  the  simulated  concentrations  of 
perchlorate  downgradient  of  the  treatment  system.  Extending  the  operating  time 
from  300  to  600  days  resulted  in  more  nondominated  solutions  achieving 
maximum  downgradient  concentrations  in  the  microgram-per-liter  range. 

•  The  technology’s  ability  to  remove  perchlorate  mass  is  poorly  correlated  with  its 
ability  to  achieve  diminished  downgradient  perchlorate  concentrations,  especially 
for  a  contaminated  site  with  relatively  low  hydraulic  conductivity  and  regional 
gradient.  In  other  words,  mass  removal  and  diminished  downgradient 
concentration  are  not  redundant  objectives.  Therefore,  decision  makers  must 
separately  consider  and  weight  each  of  these  remediation  goals  when  deciding  on 
design  parameters.  It  appears  important  to  include  downgradient  perchlorate 
concentration  as  either  an  additional  objective  or  a  constraint  when  implementing 
a  multi-objective  optimization  scheme. 

•  The  ability  of  this  technology  to  remediate  perchlorate-contaminated  groundwater 
is  very  difficult  to  ascertain  for  a  variety  of  reasons.  Results  from  the  four 
simulations  are  based  upon  numerous  simplifying  assumptions  about  kinetic 
parameters  and  site  properties,  which  may  not  be  representative  of  field 
conditions.  As  Parr  (2002)  pointed  out,  literature  values  of  kinetic  parameters  are 
“highly  variable  and  sparse.”  The  model  requires  validation  with  real-world  data 
to  properly  judge  the  capability  of  in  situ  bioremediation  technology. 
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5.3  RECOMMENDATIONS 


Take  advantage  of  the  MOGA’s  power  by  formulating  more  general  multi-objective 
problems.  Due  to  time  constraints,  this  study  placed  somewhat  artificial  (albeit  realistic) 
limits  on  several  of  the  design  parameters  that  were  being  optimized.  For  instance,  the 
study  restricted  the  number  of  treatment  wells  to  two,  and  also  imposed  restrictions  on 
well  location  (e.g.  the  wells  were  spaced  symmetrically  about  the  site  centerline  and  were 
at  a  specified  distance  downgradient  of  the  source)  and  pumping  rates  (the  two  wells  each 
pumped  at  the  same  rate).  More  generally,  HK-MOGA  could  be  applied  where  these 
decision  variables  could  all  be  optimized. 

To  avoid  over-constraining  the  problem,  CIO4'  concentration  downgradient  of  the 
treatment  system  was  not  specified  as  a  constraint.  In  reality,  downgradient 
concentration  would  be  an  important  parameter  to  use  in  defining  technology 
effectiveness.  Additionally,  results  of  this  study  confirm  that  1)  reduction  of  contaminant 
mass  and  2)  reduction  of  downgradient  contaminant  concentration  would  not  be 
redundant  objectives  for  the  in  situ  bioremediation  technology.  Future  studies  could 
implement  downgradient  concentration  as  a  constraint  or  even  as  an  additional  objective 
function  to  be  minimized. 

Accuracy  of  model  output  could  improve  by  modifying  the  perchlorate- 
contaminated  site  model  (Figure  3).  More  numerous,  smaller  cells  could  improve 
accuracy  of  the  technology  model’s  spatially-dependent  attributes  like  perchlorate  mass 
removal  and  downgradient  concentration,  although  this  modification  would  be  more 
computationally  expensive. 
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Improve  search  performance  by  experimenting  with  the  HK-MOGA’s  parameters 
or  adopting  a  different  MOGA.  There  was  no  theoretical  basis  for  selecting  values  for 
|Pop0|,  \MP\,  N,  ashare,  and  CTsharej  factors  such  as  off-line  experimentation  with  HK- 
MOGA,  literature  values,  and  resource  constraints  guided  the  author  in  assigning  values 
to  these  parameters.  Additional  simulations  could  provide  an  empirical  means  to 
determine  the  optimal  blend  of  parameter  values  that  improve  the  search  for  Pareto 
optimal  solutions. 

Chapter  2  described  several  MOGAs,  such  as  NPGA,  NSGA,  NSGA-II,  and 
SPEA.  For  this  study,  NPGA  was  selected  as  the  most  appropriate  due  to  its  ease  of 
implementation  and  its  success  in  solving  a  multi-objective  groundwater  remediation 
problem  with  objectives  similar  to  the  ones  in  this  study  (Erickson  et  al.,  2001).  Newer 
MOGAs,  like  NSGA-II,  and  SPEA,  have  an  advantage  over  NPGA  and  NSGA  because 
their  fitness  sharing  methods  do  not  require  user-specified  parameters.  However,  as  can 
be  seen  from  the  literature,  research  in  this  area  is  advancing  rapidly,  and  in  a  future 
study  a  more  efficient  algorithm  may  be  available  for  application  to  this  problem. 

Validate  the  model  with  field  data.  A  field  evaluation  of  in  situ  bioremediation  of  TCE 
using  HFTWs  has  already  been  performed  at  the  Edwards  AFB,  Site  19  (McCarty  et  al., 
1998).  A  similar  field  evaluation  is  scheduled  to  begin  this  year  to  study  a  pilot-scale 
system  to  treat  perchlorate-contaminated  groundwater.  Results  and  analysis  of  this 
evaluation  will  offer  the  opportunity  to  ascertain  kinetic  parameters  for  use  in  technology 
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design.  The  results  obtained  from  operation  of  the  pilot-scale,  in  situ  bioremediation 
system  could  also  provide  data  necessary  to  validate  Parr’s  (2002)  technology  model. 
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APPENDIX  A:  TECHNOLOGY  MODEL  EQUATIONS 


The  Parr  (2002)  technology  model  is  a  set  of  partial  differential  equations 
representing  flow  and  transport  (equations  A.l  -  A. 4),  biological  reactions  (equations  A. 5 
-  A. 8),  and  biomass  growth  (equation  A.  12)  in  a  subsurface  system  with  microorganisms 
utilizing  an  electron  donor  (acetate)  to  reduce  three  electron  acceptors  (CIO4",  oxygen, 
and  nitrate).  This  section  provides  further  detail  on  the  equations  that  are  the  technology 
model. 

FLOW  AND  TRANSPORT  MODEL  (PARR,  2002) 

The  Parr  (2002)  technology  model  solves  four  separate,  3 -dimensional 
advection/dispersion  equations  (A.l  -  A.4)  that  represent  transport  of  the  electron  donor 
(acetate  CH3COO')  and  three  electron  acceptors  (oxygen,  nitrate,  and  CIO4’), 
respectively.  The  left  side  of  equation  A.  1  includes  a  retardation  factor  ( R )  that  accounts 
for  sorption  of  the  electron  donor;  sorption  is  assumed  to  be  an  equilibrium  process  that  is 
both  linear  and  reversible.  On  the  other  hand,  equations  A.2  -  A.4  use  R=l,  as  the 
electron  acceptors  (O2,  NO3",  CIO4")  are  assumed  to  be  non-sorbing.  The  right-hand  sides 
of  equations  A.  1-A.4  include  dispersion  terms  (DV  C)  and  advection  terms  ( vV C)  for  the 
electron  donor  and  acceptors.  The  rightmost  terms  in  equations  A.l  -  A.4  are  source/sink 
tenns  that  represent  production/consumption  rates  (r)  due  to  microbial  redox  reactions. 
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where 


dc 


don 


dt 


■R=D-V2Cdon-vVCdon+r 


donor 


(A.l) 


proxy 

- - =  D-V2Coxy’-v  VCoxy  +r 

dt  oxy 


(A.2) 


dC' 


n  V7  2  /'■•f  nit  s~i  nitn  . 

=  D-v  C  -v-VC  +r, 


dt 


(A. 3) 


dC 


per 


dt 


=  D-V2Cper-v-VCper+r 


(A.4) 


O’0'1  =  concentration  of  the  electron  donor  (acetate)  (mg/L) 

C°xy  =  concentration  of  oxygen  (an  electron  acceptor)  (mg/L) 

Cmt  =  concentration  of  nitrate  (an  electron  acceptor)  (mg/L) 

Cpei  =  concentration  of  CIO4'  (an  electron  acceptor)  (mg/L) 
t  =  time  (days). 

D  =  dispersion  (mVday); 

v  =  average  linear  velocity  of  groundwater  (m/day) 

R  =  retardation  factor;  accounts  for  sorption  of  the  electron  donor 
r donor  =  rate  of  electron  donor  consumption  (mg  donor/L/day) 
roxy  =  rate  of  oxygen  consumption  (mg  oxygen/L/day) 
rmt  =  rate  of  nitrate  consumption  (mg  nitrate/L/day) 
rper  =  rate  of  CIO4"  consumption  (mg  CKV/L/day) 
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BIOLOGICAL  TREATMENT  SUBMODEL  (PARR,  2002) 


Equations  A. 5  -  A. 8  and  A.  12  comprise  the  biological  treatment  submodel.  As 
mentioned  previously,  the  rightmost  terms  in  equations  A.  1  -  A.4  are 
production/consumption  rates  due  to  microbial  redox  reactions;  these  tenns  are 
formulated  based  on  a  dual-Monod,  multi-electron  acceptor  biodegradation  model 
proposed  by  Envirogen  (Parr,  2002).  These  reaction  terms  are  defined  as  follows  (note 
that  the  negative  sign  indicates  consumption): 


dCd 


dt 


-=~X-(rc 


■hr  r  i 

don ,  oxy  don ,  nit  don ,  per  ) 


(A. 5) 


oxy 


dCoxy 
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X-Foxy-K 
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oxy 

don  _|_  don 
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+ 
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■s 

i _ 

j^oxy  _j_  oxy 

where 


Km 


Kmt  +  C 


(A. 8) 


r donor  =  rate  of  electron  donor  consumption  (mg  donor/L/day) 
r don, oxy  =  specific  rate  of  electron  donor  consumption  using  oxygen  as  an 
electron  acceptor  (mg  donor/mg  biomass/day);  see  equation  A. 9 
r don, nit  =  specific  rate  of  electron  donor  consumption  using  nitrate  as  an 
electron  acceptor  (mg  donor/mg  biomass/day);  see  equation  A.  10 
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r don, per  =  specific  rate  of  electron  donor  consumption  using  CIO4’  as  an 
electron  acceptor  (mg  donor/mg  biomass/day);  see  equation  A.  1 1 
roxy  =  rate  of  oxygen  consumption  (mg  oxygen/L/day) 
rmt  =  rate  of  nitrate  consumption  (mg  nitrate/L/day) 
rper  =  rate  of  C104‘  consumption  (mg  ClCV/L/day) 

X  =  concentration  of  active  microorganisms  (mg/L) 

Foxy  =  stoichiometric  coefficient  for  the  donor  (acetate)-oxygen  reaction 
(mg  oxygen/mg  donor)  where  the  stoichiometric  coefficient 
accounts  for  the  electron  acceptor  requirement  for  biomass 
production  based  on  the  following  stoichiometry  (C5H9NO3 
represents  the  chemical  formula  for  biomass)  (Envirogen,  2002a): 

02  +  0.64  CH3COOH  +  0.056  NH4OH  ->  0.056  C5H9N03  +  0.168  H20  +  1.0  H2C03 

Fnjt  =  stoichiometric  coefficient  for  the  donor  (acetate)-nitrate  reaction 
(mg  nitrate/mg  donor)  where  the  coefficient  accounts  for  the 
electron  acceptor  requirement  for  biomass  production  (Envirogen, 
2002a): 

N03-+  0.786  CH3COOH  <->  0.056  C5H9N03  +  0.472N2+  0.528  H20  +  0.292  H2C03  +  HC03' 

Fper  =  stoichiometric  coefficient  for  the  donor  (acetateVCICV  reaction 
(mg  CKV/mg  donor)  where  the  coefficient  accounts  for  the 
electron  acceptor  requirement  for  biomass  production  (Envirogen, 
2002a): 

CIO4  +  1.141  CH3COOH  +  0.056NH4OH  0.056  C5H9N03  +  Cl-  +  2.002  H2C03  +  0.164  H20 
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kmax  =  maximum  specific  rate  of  substrate  utilization  (mg  donor/mg 
biomass/day) 

Ks°xy  =  half  saturation  concentration  when  oxygen  (an  electron  acceptor) 
concentration  is  varied  and  limiting  (mg/L) 

Ksmt  =  half  saturation  concentration  when  nitrate  (an  electron  acceptor) 
concentration  is  varied  and  limiting  (mg/L) 

Ksl>e'  =  half  saturation  concentration  when  CIO4'  (an  electron  acceptor) 
concentration  is  varied  and  limiting  (mg/L) 

K'i'xy  =  oxygen  inhibition  coefficient  (mg/L) 

K"u  =  nitrate  inhibition  coefficient  (mg/L) 

Ksdon  =  donor  half  saturation  concentration  (mg  donor/L) 

Equations  A. 7  and  A. 8  have  inhibition  coefficients  (K,)  to  account  for  competitive 
effects  among  electron  acceptors.  Equation  A. 7  includes  an  oxygen  inhibition  coefficient 
(K,oxy)  because  the  presence  of  oxygen  inhibits  microbial  reduction  of  nitrate.  If  oxygen 
is  absent  ( C°xy  =  0),  then  K°xy  has  no  influence  on  rnit  because  the  rightmost  fraction  of 
equation  A. 7  becomes  1.  Similarly,  equation  A. 8  includes  both  oxygen  and  nitrate 
inhibition  coefficients  (K°xy  and  K,mt)  because  the  presence  of  either  oxygen  or  nitrate 
inhibits  microbial  reduction  of  ClOfi.  Parr  (2002)  cites  laboratory  results  that 
demonstrate  how  oxygen  and  nitrate  inhibit  microbial  reduction  of  less  preferable 
electron  acceptors,  like  CIO4’. 

Parr  (2002)  assumed  that  inhibition  coefficients  are  equal  to  their  respective  half 
saturation  concentrations  (i.e.  Ksoxy  =  Ki°xy  and  Ksmt  =  K.,nit).  Also,  equation  A. 5 
contains  the  tenns  rdonoxy,  rdon  mt,  and  rdonper,  which  are  defined  as  follows: 
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(A.  11) 


where 

kmax°"/oxy  =  maximum  specific  rate  of  substrate  utilization  in  the  presence 
of  oxygen  when  donor  concentration  is  varied  and  limiting  (mg 
donor/mg  biomass/day) 

kmax°n/mt  =  maximum  growth  rate  of  substrate  utilization  in  the  presence  of 
nitrate  when  donor  concentration  is  varied  and  limiting  (mg 
donor/mg  biomass/day) 

kmaxdon/per  =  maximum  specific  rate  of  substrate  utilization  in  the  presence 
of  ClOfi  when  donor  concentration  is  varied  and  limiting  (mg 
donor/mg  biomass/day) 

gdon/oxy  _  p, a )  f  saturation  concentration  of  the  electron  donor  in  the 

presence  of  oxygen  when  donor  (acetate)  concentration  is  varied 
and  limiting  (mg  donor/L) 

gdon/mt  _  pa]p  saturation  concentration  of  the  electron  donor  in  the 

presence  of  nitrate  when  donor  (acetate)  concentration  is  varied 
and  limiting  (mg  donor/L) 
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g  don/pet  _  saturation  concentration  of  the  electron  donor  in  the 

presence  of  CIO4'  when  donor  (acetate)  concentration  is  varied  and 
limiting  (mg  donor/L) 


The  microbial  growth/decay  equation  of  the  technology  model  is 


dX 

dt 


=X- 


(A.  12) 


where 


dX 

dt 


=  0;X<X 


min 


Xmin  =  minimum  biomass  concentration  (mg/L) 

Ybiomass  =  the  biomass  yield  per  mass  of  donor  consumed  (mg  biomass/mg 
electron  donor) 

b  =  biomass  decay  rate  (day'1) 

Equation  A.  12  assumes  that  biomass  concentration  will  never  decrease  below  some 
minimum  (Xmin). 
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APPENDIX  B:  DERIVATION  OF  ENERGY  CONSTANT 


In  this  study,  we  assume  that  each  HFTW  uses  a  static  mixer  to  blend  the  substrate 
(acetate)  with  water.  It  is  further  assumed  that  pumping  costs  depend  primarily  on  the 
energy  necessary  to  overcome  head  losses  from  the  static  mixer. 

We  assume  that  the  static  mixers  will  be  installed  in  very  smooth  pipes  of  diameter  D  = 
1.25  inches  (0.0318  m)  and  length  L  =  4.0  m,  which  is  the  distance  between  the  upper  and 
lower  screens  of  each  HFTW.  Also,  we  assume  that  water  has  a  kinematic  viscosity  v  = 
1.141  m  / sec  @  15°  C.  Dividing  flow  rate  Q  by  the  pipe’s  cross-sectional  area  A  =  nD~/4 
yields  velocity  V.  The  Reynolds  number  can  be  calculated  as  follows: 

Re  =  VD/v  (B.l) 

Recall  from  Section  3.6  that  the  minimum  and  maximum  values  for  Q  are  10  and  150 
m  /day,  which  correspond  to  Reynolds  numbers  of  4,060.  and  60,900,  respectively. 
Hence,  flow  is  clearly  turbulent  (Re  >  2,000)  for  all  possible  values  of  Q. 

The  Blasius  equation  can  be  used  to  obtain  the  friction  factor / for  3,000  <  Re  <  100,000 
(Daugherty  &  Franzini,  1965:212): 

/=  0.316/(Re°'25)  (B.2) 

Knowing/ enables  computation  of  the  empty-pipe  head  loss  Iil  via  the  Darcy- Weisbach 
equation  (Mays,  2001:419): 

hL  =  8fLQ2/(n2gD5)  (B.3) 
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Cleveland  Eastern  Mixers  supplied  the  following  equation  for  computing  the  flow 
coefficient  (Cf)  for  Re  >  1000: 

Cf=  -15.9  +  8.41  ln(Re)  (B.4) 

The  flow  coefficient  is  used  to  account  for  the  headloss  through  the  static  mixer  elements 
{hL, static),  using  equation  B.5  below: 

hL,static  =  Cf-hL  (B.5) 

Finally,  given  that  the  unit  weight  of  water  is  y  [units  of  force-per-m  ],  the  cnergy-per-m 
of  water  required  to  overcome  Hl, static  can  be  written  as 

E  —  y  llL, static  (B.6) 
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